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Preface 


I  have  had  the  great  privilege  of  working  in  the  time  and  frequency  field  over  the  span  of  my  career.  I  have  seen 
atomic  frequency  standards  shrink  from  racks  of  equipment  to  chip  scale,  and  be  manufactured  by  the  tens  of 
thousands,  while  primary  standards  and  the  time  dissemination  networks  that  support  them  have  improved  by  several 
orders  of  magnitude.  During  the  same  period,  significant  advances  have  been  made  in  our  ability  to  measure  and 
analyze  the  performance  of  those  devices.  This  Handbook  summarizes  the  techniques  of  frequency  stability  analysis, 
bringing  together  material  that  I  hope  will  be  useful  to  the  scientists  and  engineers  working  in  this  field. 


I  acknowledge  the  contributions  of  many  colleagues  in  the  Time  and  Frequency  community  who  have  contributed  the 
analytical  tools  that  are  so  vital  to  this  field.  In  particular,  I  wish  to  recognize  the  seminal  work  of  J. A.  Barnes  and 
D.W.  Allan  in  establishing  the  fundamentals  at  NBS,  and  D.A.  Howe  in  carrying  on  that  tradition  today  at  NIST. 
Together  with  such  people  as  M.A.  Weiss  and  C.A.  Greenhall,  the  techniques  of  frequency  stability  analysis  have 
advanced  greatly  during  the  last  45  years,  supporting  the  orders-of-magnitude  progress  made  on  frequency  standards 
and  time  dissemination. 

I  especially  thank  David  Howe  and  the  other  members  of  the  NIST  Time  and  Frequency  Division  for  their  support, 
encouragement,  and  review  of  this  Handbook. 
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1  Introduction 


This  handbook  describes  practical  techniques  for  frequency  stability  analysis.  It  covers  the  definitions  of  frequency 
stability  ,  measuring  systems  and  data  formats,  pre-processing  steps,  analysis  tools  and  methods,  post-processing  steps, 
and  reporting  suggestions.  Examples  are  included  for  many  of  these  techniques.  Some  of  the  examples  use  the 
Stable32  program  [1],  which  is  a  tool  for  studying  and  performing  frequency  stability  analyses.  Two  general 
references  [2,3]  for  this  subject  are  also  given. 

This  handbook  can  be  used  both  as  a  tutorial  and  as  a  reference.  If  this  is  your  first  exposure  to  this  field,  you  may 
find  it  helpful  to  scan  the  sections  to  gain  some  perspective  regarding  frequency  stability  analysis.  I  strongly 
recommend  consulting  the  references  as  part  of  your  study  of  this  subject  matter.  The  emphasis  is  on  time  domain 
stability  analysis,  where  specialized  statistical  variances  have  been  developed  to  characterize  clock  noise  as  a  function 
of  averaging  time.  I  present  methods  to  perform  those  calculations,  identify  noise  types,  and  determine  confidence 
limits.  It  is  often  important  to  separate  deterministic  factors  such  as  aging  and  environmental  sensitivity  from  the 
stochastic  noise  processes.  One  must  always  be  aware  of  the  possibility  of  outliers  and  other  measurement  problems 
that  can  contaminate  the  data. 

Suggested  analysis  procedures  are  recommended  to  gather  data,  preprocess  it,  analyze  stability,  and  report  results. 
Throughout  these  analyses,  it  is  worthwhile  to  remember  R.W.  Hamming's  axiom  that  "'the  purpose  of  computing  is 
insight,  not  numbers."  Analysts  should  feel  free  to  use  their  intuition  and  experiment  with  different  methods  that  can 
provide  a  deeper  understanding. 


References  for  Introduction 

1.  The  Stable32  Program  for  Frequency  Stability  Analysis,  Hamilton  Technical  Services,  Beaufort,  SC 
29907,  http://www.wriley.com. 

2.  D.B  Sullivan,  D.W  Allan,  D.A  Howe,  and  F.L  Walls,  eds.,  "Characterization  of  clocks  and  oscillators," 
Natl.  Inst.  Stand.  Technol.  Technical  Note  1337, 
http://tf.nist.gov/timefreq/general/pdf/868.pdf  (March  1990). 

3.  D.A.  Howe,  D.W.  Allan,  and  J. A.  Barnes,  "Properties  of  signal  sources  and  measurement  methods,"  Proc. 
35th  Freq.  Cont.  Symp.,  pp.  1-47, 

http://tfnist.gov/timefreq/general/pdf/554.pdf  (May  1981). 
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2      Frequency  Stability  Analysis 


The  time  domain  stability  analysis  of  a  frequency  source  is  concerned  with  characterizing  the  variables  x(t)  and  y(t) 
the  phase  (expressed  in  units  of  time  error)  and  the  fractional  frequency,  respectively.  It 
is  accomplished  with  an  array  of  phase  arrays  x,  and  frequency  data  arrays  y„  where  the 
index  i  refers  to  data  points  equally  spaced  in  time.  The  x,  values  have  units  of  time  in 
seconds,  and  the  y,  values  are  (dimensionless)  fractional  frequency,  Af/f  The  x(t)  time 
fluctuations  are  related  to  the  phase  fluctuations  by  <^  (t)  =  x(t)-27i:vo,  where  vq  is  the 
nominal  carrier  frequency  in  hertz.  Both  are  commonly  called  "phase"  to  distinguish 
them  from  the  independent  time  variable,  t.  The  data  sampling  or  measurement 


interval,  Tq,  has  units  of  seconds.  The  analysis  interval  or  period,  loosely  called 
''averaging  time",  t,  may  be  a  multiple  of  Xg  (t  =  mx^,  where  m  is  the  averaging  factor). 


The  objective  of  a 
frequency  stability  analysis 
is  to  characterize  the  phase 
and  frequency  fluctuations 
of  a  frequency  source  in 
the  time  and  frequency 
domains. 


The  goal  of  a  time  domain  stability  analysis  is  a  concise,  yet  complete,  quantitative  and  standardized  description  of 
the  phase  and  frequency  of  the  source,  including  their  nominal  values,  the  fluctuations  of  those  values,  and  their 
dependence  on  time  and  environmental  conditions. 

A  frequency  stability  analysis  is  normally  performed  on  a  single  device,  not  on  a  population  of  devices.  The  output  of 
the  device  is  generally  assumed  to  exist  indefinitely  before  and  after  the  particular  data  set  was  measured,  which  is  the 
(finite)  population  under  analysis.  A  stability  analysis  may  be  concerned  with  both  the  stochastic  (noise)  and 
deterministic  (systematic)  properties  of  the  device  under  test.  It  is  also  generally  assumed  that  the  stochastic 
characteristics  of  the  device  are  constant  (both  stationary  over  time  and  ergodic  over  their  population).  The  analysis 
may  show  that  this  is  not  true,  in  which  case  the  data  record  may  have  to  be  partitioned  to  obtain  meaningful  results.  It 
is  often  best  to  characterize  and  remove  deterministic  factors  (e.g.,  frequency  drift  and  temperature  sensitivity)  before 
analyzing  the  noise.  Environmental  effects  are  often  best  handled  by  eliminating  them  from  the  test  conditions.  It  is 
also  assumed  that  the  frequency  reference  instability  and  instrumental  effects  are  either  negligible  or  removed  from 
the  data.  A  common  problem  for  time  domain  frequency  stability  analysis  is  to  produce  results  at  the  longest  possible 
analysis  interval  in  order  to  minimize  test  time  and  cost.  Computation  time  is  generally  not  as  much  a  factor. 

2.1  Background 

The  field  of  modem  frequency  stability  analysis  began  in  the  mid  1960's  with  the  emergence  of  improved  analytical 
and  measurement  techniques.  In  particular,  new  statistics  became  available  that  were  better  suited  for  common  clock 
noises  than  the  classic  N-sample  variance,  and  better  methods  were  developed  for  high  resolution  measurements  (e.g., 
heterodyne  period  measurements  with  electronic  counters,  and  low  noise  phase  noise  measurements  with  double- 
balanced  diode  mixers).  A  seminal  conference  on  short-term  stability  in  1964  [1],  and  the  introduction  of  the  two- 
sample  (Allan)  variance  in  1966  [2]  marked  the  beginning  of  this  new  era,  which  was  summarized  in  a  special  issue  of 
the  Proceedings  of  the  IEEE  in  1966  [3].  This  period  also  marked  the  introduction  of  commercial  atomic  frequency 
standards,  increased  emphasis  on  low  phase  noise,  and  the  use  of  the  LORAN  radio  navigation  system  for  global 
precise  time  and  frequency  transfer.  The  subsequent  advances  in  the  performance  of  frequency  sources  depended 
largely  on  the  improved  ability  to  measure  and  analyze  their  stability.  These  advances  also  mean  that  the  field  of 
frequency  stability  analysis  has  become  more  complex.  It  is  the  goal  of  this  handbook  to  help  the  analyst  deal  with 
this  complexity. 

An  example  of  the  progress  that  has  been  made  in  frequency  stability  analysis  from  the  original  Allan  variance  in 
1966  through  Theol  in  2003  is  shown  in  the  plots  below.  The  error  bars  show  the  improvement  in  statistical 
confidence  for  the  same  data  set,  while  the  extension  to  longer  averaging  time  provides  better  long-term  clock 
characterization  without  the  time  and  expense  of  a  longer  data  record. 
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Figure  1.  Progress  in  frequency  stability  analysis,  (a)  Original  Allan,  (b)  Overlapping  Allan,  (c)  Total,  (d)  Theol.  (e) 
Overlapping  and  Theol. 
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This  handbook  includes  detailed  information  about  these  (and  other)  stability  measures. 


References  for  Frequency  Stability  Analysis 

1.  Proc.  of  the  IEEE-NASA  Symposium  on  the  Definition  and  Measurement  of  Short-Term  Frequency 
Stability,  NASA  SP-80,  (Nov.  1964). 

2.  D.W.  Allan,  "The  Statistics  of  Atomic  Frequency  Standards,"  Proc.  IEEE,  54(2):  221-230(Feb.  1966). 

3.  Special  Issue  on  Frequency  Stability,  Proc.  IEEE,  54(2)(Feb.  1966). 


4 


3      Definitions  and  Terminology 


Specialized  definitions  and 
terminology  are  used  for 
frequency  stability 
analysis.  


The  field  of  frequency  stability  analysis,  like  most  others,  has  its  own  specialized 
definitions  and  terminology.  The  basis  of  a  time  domain  stability  analysis  is  an  array 
of  equally  spaced  phase  (really  time  error)  or  fractional  frequency  deviation  data 
arrays,  x,  and  v„  respectively,  where  the  index  /  refers  to  data  points  in  time.  These 
data  are  equivalent,  and  conversions  between  them  are  possible.   The  x  values  have 

units  of  time  in  seconds,  and  the  v  values  are  (dimensionless)  fractional  frequency,  Af/f.  The  x(t)  time  fluctuations  are 
related  to  the  phase  fluctuations  by  4>(t)  =  x(t)  ■  Itivo,  where  i',;is  the  carrier  frequency  in  hertz.  Both  are  commonly 
called  "phase"  to  distinguish  them  from  the  independent  time  variable,  /.  The  data  sampling  or  measurement  interval. 
To,  has  units  of  seconds.  The  analysis  or  averaging  time,  r,  may  be  a  multiple  of  r,;(r  =  wr^,  where  m  is  the  averaging 
factor).  Phase  noise  is  fiindamental  to  a  frequency  stability  analysis,  and  the  type  and  magnitude  of  the  noise,  along 
with  other  factors  such  as  aging  and  environmental  sensitivity,  determine  the  stability  of  the  frequency  source. 

3.1.  Noise  Model 

A  frequency  source  has  a  sine  wave  output  signal  given  by  [1] 

V{t)  =  [V,  +  £(t)]  sm[27rv,t  +  m] ,  ( D 

where   Vq  =  nominal  peak  output  voltage 
£ftj  =  amplitude  deviation 
Vq  =  nominal  frequency 
(pOj  =  phase  deviation. 

For  the  analysis  of  frequency  stability,  we  are  concerned  primarily  with  the  ^(t)  term.  The  instantaneous  frequency  is 
the  derivative  of  the  total  phase: 

KO  =  Vo+  — (2) 
ZTT  at 

For  precision  oscillators,  we  define  the  fractional  frequency  as 

,  ,    Af    v{t)-Vr,        1    d(/>  dx 

/  Vq         IttVq  dt  dt 

where 

x{t)^m^27rv,.  (4) 

3.2.  Power  Law  Noise 

It  has  been  found  that  the  instability  of  most  frequency  sources  can  be  modeled  by  a  combination  of  power-law  noises 
having  a  spectral  density  of  their  fractional  frequency  fluctuations  of  the  form  Sy(/)  ocf,  where  /is  the  Fourier  or 
sideband  frequency  in  hertz,  and  a  is  the  power  law  exponent. 

Noise  Type  a 

White  PM  (W  PM)  2 

Flicker  PM  (F  PM)  1 

White  FM  (W  FM)  0 

Flicker  FM  (F  FM)  -1 

Random  Walk  FM  (RW  FM)  -2 
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Flicker  Walk  FM  (FW  FM) 
Random  Run  FM  (RR  FM) 


-3 
-4 


Examples  of  the  four  most  common  of  these  noises  are  shown  in  Table  1. 


Table  1.  Examples  of  the  four  most  common  noise  types. 

POWER  LAW  NOISE  SPECTRA 


WHITE  N0I9£  fs 


FLICKER  NOISE  t" 


HAWKJM  WALK  HCBSE  1-2 


FUCKER  WALK  NOISE  H 


1  I 


!  1- 


3.3.  Stability  Measures 

The  standard  measures  for  frequency  stability  in  the  time  and  frequency  domains  are  the  overlapped  Allan  deviation, 
ay(i),  and  the  SSB  phase  noise,  £(f),  as  described  in  more  detail  later  in  this  handbook. 

3.4.  Differenced  and  Integrated  Noise 

Taking  the  differences  between  adjacent  data  points  plays  an  important  role  in  frequency  stability  analysis  for 
performing  phase  to  frequency  data  conversion,  calculating  Allan  (and  related)  variances,  and  doing  noise 
identification  using  the  lag  1  autocorrelation  method  [2].  Phase  data  .Yf/)  may  be  converted  to  fractional  frequency 
data.  y(t)  by  taking  the  first  differences  x,.  /-  x,  of  the  phase  data  and  dividing  by  the  sampling  interval  r.  The  Allan 
variance  is  based  on  the  first  differences  y,./  -  y,  of  the  fractional  frequency  data  or,  equivalently,  the  second 
differences  y,^]-  2y,^  i  +  y,  of  the  phase  data.  Similarly,  the  Hadamard  variance  is  based  on  third  differences  Xr^s  - 
3x,+2  +  3x,.  i-x,  of  the  phase  data. 

Taking  the  first  differences  of  a  data  set  has  the  effect  of  making  it  less  divergent.  In  terms  of  its  spectral  density,  the 
a  value  is  increased  by  2.  For  example,  flicker  FM  data  (a  =  -1)  is  changed  into  flicker  PM  data  (a  =  +1).  That  is 
the  reason  that  the  Hadamard  variance  is  able  to  handle  more  divergent  noise  types  (a  >  -4)  than  the  Allan  variance 
(a  >  -2)  can.  It  is  also  the  basis  of  the  lag  1  autocorrelation  noise  identification  method  whereby  first  differences  are 
taken  until  a  becomes  >0.5.  The  plots  below  show  random  run  noise  differenced  first  to  random  walk  noise  and  again 
to  white  noise. 
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Random  Rui 


Random  VUak  Noise 


(a)  Original  random  run  (RR) 
noise 


1 1  1 

L  1 

400  600 


(b)  Differenced  RR  noise  = 
random  walk  (RW)  noise 


BOO  1000 


(c)  Differenced  RW  noise  = 
white  (W)  noise 

Figure  2.  (a)  Random  run  noise,  difference  to  (b)  random  walk  noise 
and  (c)  white  noise. 


The  more  divergent  noise  types  are  sometimes  referred  to  by  their  color.  White  noise  has  a  flat  spectral  density  (by 
analogy  to  white  light).  Flicker  noise  has  an  f '  spectral  density,  and  is  called  pink  or  red  (more  energy  toward  lower 
frequencies).  Continuing  the  analogy,  f"  (random  walk)  noise  is  called  brown,  and  f'  (flicker  walk)  noise  is  called 
black,  although  that  terminology  is  seldom  used  in  the  field  of  frequency  stability  analysis. 

Integration  is  the  inverse  operation  of  differencing.  Numerically  integrating  frequency  data  converts  it  into  phase  data 
(with  an  arbitrary  initial  value).  Such  integration  subtracts  2  from  the  original  a  value.  For  example,  the  random  run 
data  in  Figure  2(a)  was  generated  by  simulating  random  walk  FM  data  and  converting  it  to  phase  data  by  numerical 
integration. 

3.5.  Glossary 

See  the  Glossary  chapter  at  the  end  of  this  handbook  for  brief  defmitions  of  many  of  the  important  terms  used  in  the 
field  of  frequency  stability  analysis. 


References  for  Definitions  and  Terminology 

1.  "IEEE    Standard    Defmitions    of    Physical    Quantities    for    Fundamental    Frequency    and  Time 
Metrology-Random  instabilities,''  IEEE  Std.  11 39  (July  1999). 

2.  W.J.  Riley  and  C.A.  Greenhall,  "Power  law  noise  identification  using  the  lag  1  autocorrelation,"  Proc.  18th 
European  Frequency  and  Time  Forum,  University  of  Surrey,  Guildford,  U.K.  (April  5-7,  2004). 
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4  Standards 


Standards  have  been  adopted  for  the  measurement  and  characterization  of  frequency  stability, 
as  shown  in  the  references  below  [1-5].  These  standards  define  terminology,  measurement 
methods,  means  for  characterization  and  specification,  etc.  In  particular,  IEEE-Std-1 1 39 
contains  definitions,  recommendations,  and  examples  for  the  characterization  of  frequency 
stability. 

References  for  Standards 

1.  "Characterization  of  frequency  and  phase  noise,  Intl.  Consult.  Comm.  (C.C.I.R.),  Report  580,"  pp.  142-150 
(1986). 

2.  MIL-PRF-553 1 0,  "Oscillator,  crystal  controlled,  general  specification  for  (2006)." 

3.  R.L.  Sydnor,  ed.,  "The  selection  and  use  of  precise  frequency  systems,"  ITU-R  Handbook  (1995). 

4.  "Guide  to  the  expression  of  uncertainty  in  measurement,"  Intl.  Stand.  Org.  (ISO),  ISBN  92-67-10188-9 
(1995). 

5.  "IEEE  standard  definitions  of  physical  quantities  for  fundamental  frequency  and  time  metrology-Random 
instabilities,"  IEEE  Std.  1 139  (July  1999). 


Several  standards 
apply  to  the  field  of 
frequency  stability 
analysis.  
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5      Time  Domain  Stability 


The  stability  of  a  frequency  source  in  the  time  domain  is  based  on  the  statistics  of  its 
phase  or  frequency  fluctuations  as  a  function  of  time,  a  form  of  time  series  analysis  [1]. 
This  analysis  generally  uses  some  type  of  variance,  a  second  moment  measure  of  the 
fluctuations.  For  many  divergent  noise  types  commonly  associated  with  frequency 
sources,  the  standard  variance,  which  is  based  on  the  variations  around  the  average 
value,  is  not  convergent,  and  other  variances  have  been  developed  that  provide  a  better 
characterization  of  such  devices.  A  key  aspect  of  such  a  characterization  is  the  dependence  of  the  variance  on  the 
averaging  time  used  to  make  the  measurement,  which  dependence  shows  the  properties  of  the  noise. 


Time  domain  stability 
measures  are  based  on 
the  statistics  of  the  phase 
or  frequency  fluctuations 
as  a  function  of  time. 


5.1 


Sigma-Tau  Plots 


The  most  common  way  to  express  the  time  domain  stability  of  a  frequency  source  is  by  means  of  a  sigma-tau  plot  that 
shows  some  measure  of  frequency  stability  versus  the  time  over  which  the  frequency  is  averaged.  Log  sigma  versus 
log  tau  plots  show  the  dependence  of  stability  on  averaging  time,  and  show  both  the  stability  value  and  the  type  of 
noise.  The  power  law  noises  have  particular  slopes,  p,  as  shown  on  the  following  log  a  versus  log  x  plots,  and  a  and 
p  are  related  as  shown  in  the  table  below: 


Noise 

W  PM 

2 

-2 

F  PM 

1 

~-2 

W  FM 

0 

-1 

F  FM 

-1 

0 

RW  FM 

_2 

1 

The  log  a  versus  log  i  slopes  are  the  same  for  the  two  PM  noise  types,  but  are  different  on  a  Mod  sigma  plot,  which  is 
often  used  to  distinguish  between  them. 
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Sigma  Tau  Diagram 
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Mod  Sigma  Tau  Diagram 
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Figure  3.  (a)  Sigma  tau  diagram,  (b)  Mod  sigma  tau  diagram. 


5.2  Variances 

Variances  are  used  to  characterize  the  fluctuations  of  a  frequency  source  [2-3].  These  are  second-moment  measures 
of  scatter,  much  as  the  standard  variance  is  used  to  quantify  the  variations  in,  say,  the  length  of  rods  around  a  nominal 
value.  The  variations  from  the  mean  are  squared,  summed,  and  divided  by  one  less  than  the  number  of  measurements; 
this  number  is  called  the  "degrees  of  freedom." 

Several  statistical  variances  are  available  to  the  frequency  stability  analyst,  and  this  section  provides  an  overview  of 
them,  with  more  details  to  follow.  The  Allan  variance  is  the  most  common  time  domain  measure  of  frequency 
stability,  and  there  are  several  versions  of  it  that  provide  better  statistical  confidence,  can  distinguish  between  white 
and  flicker  phase  noise,  and  can  describe  time  stability.  The  Hadamard  variance  can  better  handle  frequency  drift  and 
more  divergence  noise  types,  and  several  versions  of  it  are  also  available.  The  newer  Total  and  Theol  variances  can 
provide  better  confidence  at  longer  averaging  factors. 

There  are  two  categories  of  stability  variances:  unmodified  variances,  which  use  d'*'  differences  of  phase  samples,  and 
modified  variances,  which  use  d*^  differences  of  averaged  phase  samples.  The  Allan  variances  correspond  to  d  =  2, 
and  the  Hadamard  variances  to  d  =  3.  The  corresponding  variances  are  defined  as  a  scaling  factor  times  the  expected 
value  of  the  differences  squared.  One  obtains  unbiased  estimates  of  this  variance  from  available  phase  data  by 
computing  time  averages  of  the  differences  squared.  The  usual  choices  for  the  increment  between  estimates  (the  time 
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step)  are  the  sample  period  tq  and  the  analysis  period  t,  a  multiple  of  Tq.  These  give  respectively  the  overlapped 
estimator  and  non-overlapped  estimators  of  the  stability. 


Variance  Type 

Standard 
Allan 

Overlapping  Allan 
Modified  Allan 
Time 
Hadamard 

Overlapping  Hadamard 
Total 

Modified  Total 
Time  Total 
Hadamard  Total 
Theol 
TheoH 


Characteristics 

Non-convergent  for  some  clock  noises  -  don't  use 

Classic  -  use  only  if  required  -  relatively  poor  confidence 

General  purpose  -  most  widely  used  -  first  choice 

Used  to  distinguish  W  and  F  PM 

Based  on  modified  Allan  variance 

Rejects  frequency  drift,  and  handles  divergent  noise 

Better  confidence  than  normal  Hadamard 

Better  confidence  at  long  averages  for  Allan 

Better  confidence  at  long  averages  for  modified  Allan 

Better  confidence  at  long  averages  for  time 

Better  confidence  at  long  averages  for  Hadamard 

Provides  information  over  nearly  full  record  length 

Hybrid  of  Allan  and  TheoBR  (bias-removed  Theol)  variances 


•  All  are  second  moment  measures  of  dispersion  -  scatter  or  instability  of  frequency  from  central  value. 

•  All  are  usually  expressed  as  deviations. 

•  All  are  normalized  to  standard  variance  for  white  FM  noise. 

•  All  except  standard  variance  converge  for  common  clock  noises. 

•  Modified  types  have  additional  phase  averaging  that  can  distinguish  W  and  F  PM  noises. 

•  Time  variances  based  on  modified  types. 

•  Hadamard  ty  pes  also  converge  for  FW  and  RR  FM  noise. 

•  Overlapping  ty  pes  provide  better  confidence  than  classic  Allan  variance. 

•  Total  types  provide  better  confidence  than  corresponding  overlapping  types. 

•  TheoH  (hybrid-TheoBR)  and  Theol  (Theoretical  Variance  #1)  provide  stability  data  out  to  75  %  of  record 
length. 

•  Some  are  quite  computationally  intensive,  especially  if  results  are  wanted  at  all  (or  many)  analysis  intervals 
(averaging  times),  t.  Use  octave  or  decade  z  intervals. 


The  modified  Allan  deviation  (MDEV)  can  be  used  to  distinguish  between  white  and  flicker  PM  noise.  For  example, 
the  W  and  F  PM  noise  slopes  are  both  «  1.0  on  the  Allan  Deviation  (ADEV)  plots  in  Figure  4,  but  they  can  be 
distinguished  as  -1 .5  and  -1 .0,  respectively,  on  the  MDEV  plots. 
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Figure  4.  (a)  Slope  of  W  PM  using  Adev,  (b)  slope  of  F  PM  using  ADEV,  (c)  slope  of  W  PM  using  MDEV,  and  (d)  slope  of 
F  PM  using  MDEV. 

The  Hadamard  deviation  may  be  used  to  reject  linear  frequency  drift  when  a  stability  analysis  is  performed.  For 
example,  the  simulated  frequency  data  for  a  rubidium  frequency  standard  in  Figure  5(a)  shows  significant  drift.  Allan 
deviation  plots  for  these  data  are  shown  in  Figure  5(c)  and  (d)  for  the  original  and  drift-removed  data.  Notice  that, 
without  drift  removal,  the  Allan  deviation  plot  has  a  +i  dependence  at  long  i,  a  sign  of  linear  frequency  drift. 
However,  as  seen  in  Figure  5(b),  the  Hadamard  deviation  for  the  original  data  is  nearly  the  same  as  the  Allan 
deviation  after  drift  removal,  but  it  has  lower  confidence  for  a  given  x. 
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Figure  5.  (a)  Simulated  frequency  data  for  a  rubidium  frequency  standard,  (b)  overlapping  Hadamard  with  drift,  (c) 
overlapping  sigma  with  drift,  and  (d)  overlapping  sigma  without  drift. 


References  for  Time  Domain  Stability 


1.  G.E.P.  Box  and  G.M.  Jenkins,  Time  Series  Analysis:  Forecasting  and  Control,  Holden-Day,  San  Francisco 
(1970). 

2.  J.  Rutman,  ''Characterization  of  phase  and  frequency  instabilities  in  precision  frequency  sources:  Fifteen 
years  of  progress,"  Proc.  IEEE,  66(9):  1048-1075  (1978) 

3.  S.R.  Stein.  "Frequency  and  time:  Their  measurement  and  characterization,"  Precision  Frequency  Control, 
2,  E.A.  Gerber  and  A.  Ballato,  eds.,  Academic  Press,  New  York,  ISBN  0-12-280602-6  (1985). 


The  standard  variance  should 
not  be  used  for  the  analysis  of 
frequency  stability.  


5.2.1.         Standard  Variance 

The  classic  N-sample  or  standard  variance  is  defined  as 
-—tU-^-r.  (5) 

-      1  ^ 

where  thej',  are  the  TV  fractional  frequency  values,  andj^  =  —  /  j^,  is  the  average  frequency.  The  standard  variance  is 

usually  expressed  as  its  square  root,  the  standard  deviation,  s.  It  is  not  recommended  as  a  measure  of  frequency 
stability  because  it  is  non-  convergent  for  some  types  of  noise  commonly  found  in  frequency  sources,  as  shown  in  the 
figure  below. 
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Figure  1.  Convergence  of  standard  and  Allan  deviation  for  FM  noise. 

The  standard  deviation  (upper  curve)  increases  with  the  number  of  samples  of  flicker  FM  noise  used  to  determine  it, 
while  the  Allan  deviation  (lower  curve  and  discussed  below)  is  essentially  constant. 

The  problem  with  the  standard  variance  stems  from  its  use  of  the  deviations  from  the  average,  which  is  not  stationary 
for  the  more  divergence  noise  types.  That  problem  can  be  solved  by  instead  using  the  first  differences  of  the 
fractional  frequency  values  (the  second  differences  of  the  phase),  as  described  for  the  Allan  variance  in  Section  5.2.2. 

In  the  context  of  frequency  stability  analysis,  the  standard  variance  is  used  primarily  in  the  calculation  of  the  Bl  ratio 
for  noise  recognition. 


Reference  for  Standard  Variance 

1.    D.W.  Allan,  "Should  the  Classical  Variance  be  used  as  a  Basic  Measure  in  Standards  Metrology?"  IEEE  Trans. 
Imtrum.  Meas.,  IM-36:  646-654  (1987) 


5.2.2. 


Allan  Variance 


The  Allan  variance  is  the  most  common  time  domain  measure  of  frequency  stability.  Similar  to  the  standard  variance, 
it  is  a  measure  of  the  fractional  frequency  fluctuations,  but  has  the  advantage  of  being  convergent  for  most  types  of 
clock  noise.  There  are  several  versions  of  the  Allan  variance  that  provide  better  statistical  confidence,  can  distinguish 
between  white  and  flicker  phase  noise,  and  can  describe  time  stability. 


The  original  non-overlapped  Allan,  or  two-sample  variance,  AVAR,  is  the  standard 
time  domain  measure  of  frequency  stability  [1,  2].  It  is  defined  as 


1 


M- 


2(M-l)t 


The  original  Allan  variance 
has  been  largely  superseded 
by  its  overlapping  version. 


(6) 


where  7,  is  the  /th  of  M  fractional  frequency  values  averaged  over  the  measurement  (sampling)  interval,  r.  Note  that 
these    symbols  are  sometimes  shown  with  a  bar  over  them  to  denote  the  averaging. 
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In  terms  of  phase  data,  the  Allan  variance  may  be  calculated  as 


I  /V-2 

2(N-2)t'  ,j 
where  x,  is  the  /th  of  the  N  =  M+1  phase  values  spaced  by  the  measurement  interval  r. 

The  result  is  usually  expressed  as  the  square  root,  ay(T),  the  Allan  deviation,  ADEV.  The  Allan  variance  is  the  same 
as  the  ordinar>'  variance  for  white  FM  noise,  but  has  the  advantage,  for  more  divergent  noise  types  such  as  flicker 
noise,  of  converging  to  a  value  that  is  independent  on  the  number  of  samples.  The  confidence  interval  of  an  Allan 
deviation  estimate  is  also  dependent  on  the  noise  type,  but  is  often  estimated  as  ±Oy{x)l^. 

5.2.2.         Overlapping  Samples 

Some  stabilit\  calculations  can  utilize  (fully)  overlapping  samples,  whereby  the 
calculation  is  performed  by  utilizing  all  possible  combinations  of  the  data  set,  as 
shown  in  the  diagram  and  formulae  below.  The  use  of  overlapping  samples 
improves  the  confidence  of  the  resulting  stability  estimate,  but  at  the  expense  of 
greater  computational  time.  The  overlapping  samples  are  not  completely  independent,  but  do  increase  the  effective 
number  of  degrees  of  freedom.  The  choice  of  overlapping  samples  applies  to  the  Allan  and  Hadamard  variances. 
Other  variances  (e.g.,  total)  always  use  them. 

Overlapping  samples  don't  apply  at  the  basic  measurement  interval,  which  should  be  as  short  as  practical  to  support  a 
large  number  of  overlaps  at  longer  averaging  times. 


Overlapping  samples  are  used 
to  improve  the  confidence  of 
a  stability  estimate.  


Averaging  Factor,  m  =3 
1  2 


Non-Overlapping  Samples 
3  4 


•     •  •  • 

•     •  • 


Overlapping  Samples 


Non-Overlapped  Allan        '^^^  ~  2[M -\)^^'^'*^  '^'^ 
Variance:  Stride  =  i  = 
averaging  period  = 


M-lm+\  j  +  m-\ 


2m' [M  -  2m  + 


7=1  /=7 


(8) 


(9) 


Overlapped  Allan 
Variance:  Stride  =  Tq  = 
sample  period 

Figure  7.  Comparison  of  non-overlapping  and  overlapping  sampling. 


The  following  plots  show  the  significant  reduction  in  variability,  hence  increased  statistical  confidence,  obtained  by 
using  overlapping  samples  in  the  calculation  of  the  Hadamard  deviation. 
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Figure  8.  The  reduction  in  variability  by  using  overlapping  samples  in  calculating  the  Hadamard  deviation. 


5.2.4. 


Overlapping  Allan  Variance 


The  fully  overlapping  Allan  variance,  or  AVAR,  is  a  form  of  the  normal  Allan 
variance,  cry(T),  that  makes  maximum  use  of  a  data  set  by  forming  all  possible 
overlapping  samples  at  each  averaging  time  t.  It  can  be  estimated  from  a  set  of  M 
frequency  measurements  for  averaging  time  x  =  mt,,,  where  m  is  the  averaging  factor 
and  To  is  the  basic  measurement  interval,  by  the  expression 

I  M -2m+\  f  /+m-l  1  ^ 

^li^)  =  ^  ^ — ^  Z    lL^y..r„-y^\  ■  (10) 

2m"(M-2m  +  l)   7"/    [  ^  J 


The  overlapped  Allan 
deviation  is  the  most  common 
measure  of  time-domain 
frequency  stability.  The  term 
AVAR  has  come  to  be  used 
mainly  for  this  form  of  the 
Allan  variance,  and  ADEV 
for  its  square  root.  


This  formula  is  seldom  used  for  large  data  sets  because  of  the  computationally  intensive  inner  summation.  In  terms  of 
phase  data,  the  overlapping  Allan  variance  can  be  estimated  from  a  set  of  N  =  M  +  1  time  measurements  as 


1 


2{N-2m)T 


N-2m 


(11) 


(=1 


Fractional  frequency  data,  v„  can  be  first  integrated  to  use  this  faster  formula.  The  result  is  usually  expressed  as  the 
square  root,  ay(T),  the  Allan  deviation,  ADEV.  The  confidence  interval  of  an  overlapping  Allan  deviation  estimate  is 
better  than  that  of  a  normal  Allan  variance  estimation  because,  even  though  the  additional  overlapping  differences  are 
not  all  statistically  independent,  they  nevertheless  increase  the  number  of  degrees  of  freedom  and  thus  improve  the 
confidence  in  the  estimation.  Analytical  methods  are  available  for  calculating  the  number  of  degrees  of  freedom  for 
an  estimation  of  overlapping  Allan  variance,  and  using  that  to  establish  single-  or  double-sided  confidence  intervals 
for  the  estimate  with  a  certain  confidence  factor,  based  on  Chi-squared  statistics. 


Sample  variances  are  distributed  according  to  the  expression 
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r 


df- 


(12) 


where  X'  is  the  Chi-square,  is  the  sample  variance,  cr'  is  the  true  variance,  and  df  is  the  number  of  degrees  of 
freedom  (not  necessarily  an  integer).  For  a  particular  statistic,  df'xs  determined  by  the  number  of  data  points  and  the 
noise  ty  pe. 


References  for  Allan  Variance 


D.W.  Allan,  "The  statistics  of  atomic  frequency  standards,"  Proc.  IEEE,  54(2):  221-230  (Feb.  1966). 
D.W.  Allan,  "Allan  variance, "  http://www.allanstinie.coiii  [2008]. 

"Characterization  of  frequency  stability,"  Nat.  Bur.  Stand.  (U.S.)  Tech  Note  394  (Oct.  1970). 
J. A.  Barnes,  A.R.  Chi,  L.S.  Cutler,  D.J.  Healey,  D.B.  Leeson,  T.E.  McGunigal,  J.A.  Mullen,  Jr.,  W.L.  Smith, 
R.L.  Sydnor,  R.F.C.  Vessot,  and  G.M.R.  Winkler,  "Characterization  of  frequency  stability,"  IEEE  Trans. 
Instrum.  Meas.,  20(2):  105-120  (May  1971 ) 

J.A.  Barnes,  "Variances  based  on  data  with  dead  time  between  the  measurements,"  Natl.  Inst.  Stand. 
Technol.  Technical  Note  1318  (1 990). 

C.A.  Greenhall,  "Does  Allan  variance  determine  the  spectrum?"  Proc.  1997  Intl.  Freq.  Cont.  Symp.,  pp. 
358-365  (June  1997)/ 

C.A.  Greenhall,  "Spectral  ambiguity  of  Allan  variance,"  IEEE  Trans.  Instrum.  Meas.,  47(3):  623-627  (June 
1998). 


5.2.5 


Modified  Allan  Variance 


The  modified  Allan  variance.  Mod  a-y(i).  MVAR,  is  another  common  time 
domain  measure  of  frequency  stability  [1].  It  is  estimated  from  a  set  of  M 
frequency  measurements  for  averaging  time  t  =  mto,  where  m  is  the  averaging 
factor  and  tq  is  the  basic  measurement  interval,  by  the  expression 


Use  the  modified  Allan  deviation 
to  distinguish  between  white  and 
flicker  PM  noise. 


Modd'iz) 


1 


2m  {M  -3m  +  2) 


M-2m+2  I  /+"'-! i+m~] 


,=j    V  k=i 


I 


J 


(13) 


In  terms  of  phase  data,  the  modified  Allan  variance  is  estimated  from  a  set  of  N  =  M  +  1  time  measurements  as 
1 


.V-3m+l  7+m-l 

2w-r"(iV-3m  +  l)  7-f   [  ^ 


(14) 


The  result  is  usually  expressed  as  the  square  root.  Mod  ay(x),  the  modified  Allan  deviation.  The  modified  Allan 
variance  is  the  same  as  the  normal  Allan  variance  for  m  =  1.  It  includes  an  additional  phase  averaging  operation,  and 
has  the  advantage  of  being  able  to  distinguish  between  white  and  flicker  PM  noise.  The  confidence  interval  of  a 
modified  Allan  deviation  determination  is  also  dependent  on  the  noise  type,  but  is  often  estimated  as  ±ay(T)/VN. 
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References  for  Modified  Allan  Variance 

1.  D.W.  Allan  and  J. A.  Barnes,  "A  modified  Allan  variance  with  increased  oscillator  characterization  ability," 
Proc.  35th  Freq.  Cont.  Symp.  pp.  470-474  (May  1981 ). 

2.  P.  Lesage  and  T.  Ayi,  "Characterization  of  frequency  stability:  Analysis  of  the  modified  Allan  variance  and 
properties  of  its  estimate,"  IEEE  Trans.  lustrum.  Meas.,  33(4):  332-336  (Dec.  1984). 

3.  C.A.  Greenhall,  "Estimating  the  modified  Allan  variance,"  Proc.  IEEE  1995  Freq.  Cont.  Symp.,  pp.  346-353 
(May  1995). 

4.  C.A.  Greenhall,  "The  third-difference  approach  to  modified  Allan  variance,"  IEEE  Trans.  Instrum.  Meas., 
46(3):  696-703  (June  1997). 


5.2.6. 


Time  Variance 


The  time  Allan  variance,  TVAR,  with  square  root  TDEV,  is  a  measure  of  time 
stability  based  on  the  modified  Allan  variance  [1].  It  is  defined  as 


cr;{T) 


\  J 


Mod<j]{T) 


(15) 


Use  the  time  deviation  to 
characterize  the  time  error  of 
a  time  source  (clock)  or 
distribution  system.  


In  simple  terms,  TDEV  is  MDEV  whose  slope  on  a  log-log  plot  is  transposed  by  +1  and  normalized  by  V3.  The  time 
Allan  variance  is  equal  to  the  standard  variance  of  the  time  deviations  for  white  PM  noise.  It  is  particularly  useful  for 
measuring  the  stability  of  a  time  distribution  network. 

It  can  be  convenient  to  include  TDEV  information  on  a  MDEV  plot  by  adding  lines  of  constant  TDEV,  as  shown  in 
Figure  9: 
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Figure  9.  Plot  of  MDEV  with  lines  of  constant  TDEV. 
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References  for  Time  Variance 

1.  D.W.  Allan.  D.D.  Davis.  J.  Levine.  M.A.  Weiss,  N.  Hironaka,  and  D.  Okayama.  "New  inexpensive  frequency 
calibration  service  from  Natl.  Inst.  Stand.  Technol.,"  Proc.  44th  Freq.  Cont.  Symp.,  pp.  107-1 16  (June  1990). 

2.  D.W.  Allan,  M.A.  Weiss,  and  J.L.  Jespersen.  "A  frequency-domain  view  of  time-domain  characterization  of 
clocks  and  time  and  frequency  distribution  systems,"'  Proc.  45th  Freq.  Cont.  Symp.,  pp.  667-678  (May 
1991). 


5.2.7.         Time  Error  Prediction 

The  time  error  of  a  clock  driven  b\  a  frequency  source  is  a  relatively  simple  function 
of  the  initial  time  offset,  the  frequency  offset,  and  the  subsequent  frequency  drift, 
plus  the  effect  of  noise,  as  shown  in  the  following  expression: 


AT  =  To  +  (Af/f)  ■  t  +  '/2  D  •  t-  +  cTx(t), 


(16) 


The  time  error  of  a  clock  can 
be  predicted  from  its  time  and 
frequency  offsets,  frequency 
drift,  and  noise.  


where  AT  is  the  total  time  error.  To  is  the  initial  synchronization  error,  Af/f  is  the  sum  of  the  initial  and  average 
environmentally  induced  frequencv  offsets.  D  is  the  frequency  drift  (aging  rate),  and  ax(t)  is  the  root-mean-square 
(rms)  noise-induced  time  deviation.  For  consistency,  units  of  dimensionless  fractional  frequency  and  seconds  should 
be  used  throughout. 

Because  of  the  many  factors,  conditions,  and  assumptions  involved,  and  their  variability,  clock  error  prediction  is 
seldom  easy  or  exact,  and  it  is  usually  necessary  to  generate  a  timing  error  budget. 

•    Initial  Synchronization 

The  effect  of  an  initial  time  (synchronization)  error.  To.  is  a  constant  time  offset  due  to  the  time  reference,  the  finite 
measurement  resolution,  and  measurement  noise.  The  measurement  resolution  and  noise  depends  on  the  averaging 
time. 


•     Initial  Syntonization 

The  effect  of  an  initial  frequency  (syntonization)  error,  Af/f  ,  is  a  linear  time  error.  Without  occasional 
resyntonization  (frequency  recalibration).  frequency  aging  can  cause  this  to  be  the  biggest  contributor  toward  clock 
error  for  many  frequency  sources  (e.g.,  quartz  crvstal  oscillators  and  rubidium  gas  cell  standards).  Therefore,  it  can  be 
important  to  have  a  means  for  periodic  clock  syntonization  (e.g.,  GPS  or  cesium  beam  standard).  In  that  case,  the 
syntonization  error  is  subject  to  uncertaint\  due  to  the  frequency  reference,  the  measurement  and  tuning  resolution, 
and  noise  considerations.  The  measurement  noise  can  be  estimated  by  the  square  root  of  the  sum  of  the  Allan 
variances  of  the  clock  and  reference  over  the  measurement  interval.  The  initial  syntonization  should  be  performed,  to 
the  greatest  extent  possible,  under  the  same  environmental  conditions  (e.g..  temperature)  as  expected  during 
subsequent  operation. 


•     Environmental  Sensitivity 

After  initial  syntonization,  environmental  sensitivity  is  likely  to  be  the  largest  contributor  to  time  error. 
Environmental  frequency  sensitivity  obviously  depends  on  the  properties  of  the  device  and  its  operating  conditions. 
When  performing  a  frequency  stability  analysis,  it  is  important  to  separate  the  deterministic  environmental 
sensitivities  from  the  stochastic  noise.  This  requires  a  good  understanding  of  both  the  device  and  its  environment. 
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Reference  for  Time  Error  Prediction 

D.W.  Allan  and  H.  Hellwig,  "Time  Deviation  and  Time  Prediction  Error  for  Clock  Specification,  Characterization, 
and  Application",  Proceedings  of  the  Position  Location  and  Navigation  Symposium  (PLANS),  29-36,  1978. 


5.2.8.         Hadamard  Variance 


Use  the  Hadamard  variance  to 
characterize  frequency 
sources  with  divergent  noise 
and/or  frequency  drift.  


The  Hadamard  [1]  variance  is  based  on  the  Hadamard  transform  [2],  which  was 
adapted  by  Baugh  as  the  basis  of  a  time-domain  measure  of  frequency  stability  [3]. 
As  a  spectral  estimator,  the  Hadamard  transform  has  higher  resolution  than  the  Allan 
variance,  since  the  equivalent  noise  bandwidth  of  the  Hadamard  and  Allan  spectral 
windows  are  1.2337N''t"'  and  0.476i"',  respectively  [4].  For  the  purposes  of  time- 
domain  frequency  stability  characterization,  the  most  important  advantage  of  the  Hadamard  variance  is  its 
insensitivity  to  linear  frequency  drift,  making  it  particularly  useful  for  the  analysis  of  rubidium  atomic  clocks  [5,6].  It 
has  also  been  used  as  one  of  the  components  of  a  time-domain  multivariance  analysis  [7],  and  is  related  to  the  third 
structure  function  of  phase  noise  [8]. 

Because  the  Hadamard  variance  examines  the  second  difference  of  the  fractional  frequencies  (the  third  difference  of 
the  phase  variations),  it  converges  for  the  Flicker  Walk  FM  (a  =  -3)  and  Random  Run  FM  (ot  =  -4)  power-law  noise 
types.  It  is  also  unaffected  by  linear  frequency  drift. 

For  frequency  data,  the  Hadamard  variance  is  defmed  as: 

1  M-2 

where  y,  is  the  ith  of  M  fractional  frequency  values  at  averaging  time  t. 
For  phase  data,  the  Hadamard  variance  is  defined  as: 

^'^U^'>=  r  ^  S [^-3  - ^^,.2  +  3x,,,  - Xj-  ,  (18) 

6r"(/v -3m)  TTr 
where  x,  is  the  ith  of  N  =  M  +  1  phase  values  at  averaging  time  x. 

Like  the  Allan  variance,  the  Hadamard  variance  is  usually  expressed  as  its  square-root,  the  Hadamard  deviation, 
HDEV  or  Hay(T). 


5.2.9.         Overlapping  Hadamard  Variance 

In  the  same  way  that  the  overlapping  Allan  variance  makes  maximum  use  of  a  data  The  overlapping  Hadamard 
set  by  forming  all  possible  fully  overlapping  2-sample  pairs  at  each  averaging  time  t  variance  provides  better 
the  overlapping  Hadamard  variance  uses  all  3-sample  combinations  [9].  It  can  be  confidence  than  the  non- 
estimated  from  a  set  of  M  frequency  measurements  for  averaging  time  t  =  mto  where  overlapping  version.  

m  is  the  averaging  factor  and  lois  the  basic  measurement  interval,  by  the  expression: 
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I  A/-3m  +  l  7+m-I 

^^'^^)=7-^77r";  77^  Z      Z  [>'-2.-2>',..+>^J    ,  (19) 

6m-{M -3m  +  \)T-    7~f    [  J 

where  yi  is  the  ith  of  M  fractional  frequency  values  at  each  measurement  time. 

In  terms  of  phase  data,  the  overlapping  Hadamard  variance  can  be  estimated  from  a  set  of  N  =  M  +  1  time 
measurements  as: 

Hcj]{T)^         ^      .  Z  K3,„-3x,.2™+3x_-xj\  (20) 
6(A/-j»w)r"  TTT 

where  x,  is  the  ith  of  N  =  M  +  1  phase  values  at  each  measurement  time. 

Computation  of  the  overlapping  Hadamard  variance  is  more  efficient  for  phase  data,  where  the  averaging  is 
accomplished  by  simply  choosing  the  appropriate  interval.  For  frequency  data,  an  inner  averaging  loop  over  m 
frequency  values  is  necessary.  The  resuh  is  usually  expressed  as  the  square  root.  Hay(T),  the  Hadamard  deviation, 
HDEV.  The  expected  value  of  the  overlapping  statistic  is  the  same  as  the  normal  one  described  above,  but  the 
confidence  interval  of  the  estimation  is  better.  Even  though  not  all  the  additional  overlapping  differences  are 
statistically  independent,  they  nevertheless  increase  the  number  of  degrees  of  freedom  and  thus  improve  the 
confidence  in  the  estimation.  Analytical  methods  are  available  for  calculating  the  number  of  degrees  of  freedom  for 
an  overlapping  Allan  variance  estimation,  and  that  same  theory  can  be  used  to  establish  reasonable  single-  or  double- 
sided  confidence  intervals  for  an  overlapping  Hadamard  variance  estimate  with  a  certain  confidence  factor,  based  on 
Chi-squared  statistics. 

Sample  variances  are  distributed  according  to  the  expression: 

X\p,df)  =  ^,  (21) 

where  X'  is  the  Chi-square  value  for  probability  p  and  degrees  of  freedom  df,  s-  is  the  sample  variance,  a-  is  the  true 
variance,  and  df  is  the  number  of  degrees  of  freedom  (not  necessarily  an  integer).  The  df  is  determined  by  the  number 
of  data  points  and  the  noise  tv  pe.  Given  the  df.  the  confidence  limits  around  the  measured  sample  variance  are  given 
by: 

,        {s'-df)             ,  (s'-df) 
^.,n  =    2,  and  a,,,  =  —  — .  (22) 


5.2.10.        Modified  Hadamard  Variance 

By  similarity  to  the  modified  Allan  variance,  a  modified  version  of  the  Hadamard  variance  can  be  defined  [15]  that 
employs  averaging  of  the  phase  data  over  the  m  adjacent  samples  that  define  the  analysis  r  =  m-ZQ.  In  terms  of  phase 
data,  the  three-sample  modified  Hadamard  variance  is  defined  as: 


^-4™+!  /+m-l 


2 


6m'T'[N  -Am  +  \\ 


Modc7i,{T)=    "  ^   ,  .  ^  (23) 
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where  A'^  is  the  number  of  phase  data  points  .y,  at  the  sampUng  interval  ro,  and  m  is  the  averaging  factor,  which  can 
extend  from  1  to  LiVMj.  This  is  an  unbiased  estimator  of  the  modified  Hadamard  variance,  MHVAR.  Expressions  for 
the  equivalent  number  of  x'  degrees  of  freedom  (edf)  required  to  set  MHVAR  confidence  limits  are  available  in  [2]. 


Clock  noise  (and  other  noise  processes)  can  be  described  in  terms  of  power  spectral  density,  which  can  be  modeled  as 
a  power  law  function  S  cc  where  /is  Fourier  frequency  and  a  is  the  power  law  exponent.  When  a  variance  such  as 
MHVAR  is  plotted  on  log-log  axes  versus  averaging  time,  the  various  power  law  noises  correspond  to  particular 
slopes  //.  MHVAR  was  developed  in  Reference  [15]  for  determining  the  power  law  noise  type  of  Internet  traffic 
statistics,  where  it  was  found  to  be  slightly  better  for  that  purpose  than  the  modified  Allan  variance,  MVAR,  when 
there  were  a  sufficient  number  of  data  points.  MHVAR  could  also  be  useful  for  frequency  stability  analysis,  perhaps 
in  cases  where  it  was  necessary  to  distinguish  between  short-term  white  and  flicker  PM  noise  in  the  presence  of  more 
divergent  (a=  -3  and  -4)  flicker  walk  and  random  run  FM  noises.  The  Mod  o'ni'^)  log-log  slope  /iis  related  to  the 
power  law  noise  exponent  by  /j  =  -3  -  a. 

The  modified  Hadamard  variance  concept  can  be  generalized  to  subsume  AVAR,  HVAR,  MVAR,  MHVAR,  and 
MHVARs  using  higher-order  differences: 


where  d  =  phase  differencing  order;  d  =  2  corresponds  to  MAVAR,  <:/=  3  to  MHVAR;  higher-order  differencing  is  not 
commonly  used  in  the  field  of  frequency  stability  analysis.  The  unmodified,  nonoverlapped  AVAR  and  HVAR 
variances  are  given  by  setting  m  =  1 .  The  allowable  power  law  exponent  for  convergence  of  the  variance  is  equal  to  a 
>  1  -  2d,  so  the  second  difference  Allan  variances  can  be  used  for  a  >  -3  and  the  third  difference  Hadamard 
variances  for  a  >  -5. 

Confidence  intervals  for  the  modified  Hadamard  variance  can  be  determined  by  use  of  the  edf  values  of  reference 


(24) 


dlm-T-[N  -(d  +  \)m  +  \] 


[16]. 
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5.2.11 


Total  Variance 


The  total  variance  offers 
improved  confidence  at  large 
averaging  factor  by  extending 
the  data  set  by  reflection  at 
both  ends. 


The  total  variance,  TOTVAR,  is  similar  to  the  two-sample  or  Allan  variance  and  has 
the  same  expected  value,  but  offers  improved  confidence  at  long  averaging  times  [1- 
5].  The  work  on  total  variance  began  with  the  realization  that  the  Allan  variance  can 
"collapse"  at  long  averaging  factors  because  of  symmetry  in  the  data.  An  early  idea 
was  to  shift  the  data  by  1/4  of  the  record  length  and  average  the  two  resulting  Allan 
variances.  The  next  step  was  to  wrap  the  data  in  a  circular  fashion  and  calculate  the 

average  of  all  the  Allan  variances  at  every  basic  measurement  interval,  t„.  This  technique  is  very  effective  in 
improving  the  confidence  at  long  averaging  factors  but  requires  end  matching  of  the  data.  A  further  improvement  of 
the  total  variance  concept  was  to  extend  the  data  by  reflection,  first  at  one  end  of  the  record  and  then  at  both  ends. 
This  latest  technique,  called  TOTVAR,  gives  a  very  significant  confidence  advantage  at  long  averaging  times,  exactly 
decomposes  the  classical  standard  variance  [6],  and  is  an  important  new  general  statistical  tool.  TOTVAR  is  defined 
for  phase  data  as: 

1  2 

rorvar(r)^  +^'-^^]  '  ^^^^ 
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where  r  =  mTo,  and  the  phase  values  x  measured  at  t  =  To  are  extended  by  reflection  about  both  endpoints  to  form  a 
virtual  sequence  x*  from  /  =  3-N  to  /  =  2N-2  of  length  iA^— ^.  The  original  data  are  in  the  center  of  x*  with  i  =1  to  N 
and  x*=x.  The  reflected  portions  added  at  each  end  extend  from  j  =  1  to  N-2  where  x*i_j  =  2xi-Xi+j  and  x*n+j  =  2xn- 

Xn-j- 

Totvar  can  also  be  defined  for  frequency  data  as: 

""'^^'^'^^^SC-^''-'-^-]'  '''' 

where  the  M  =  N-1  fractional  frequency  values,  y,  measured  at  t  =  Tq  (N  phase  values)  are  extended  by  reflection  at 
both  ends  to  form  a  virtual  array  y*.  The  original  data  are  in  the  center,  where  y*i  =  y,  for  i  =  1  to  M,  and  the  extended 
data  for  j  =  1  to  M-1  are  equal  to  y*i.j  =  yj  and  y*M+i  =  Ym+i-j- 

The  result  is  usually  expressed  as  the  square  root,  atotai(x),  the  total  deviation,  TOTDEV.  When  calculated  by  use  of 
the  doubly  reflected  method  described  above,  the  expected  value  of  TOTVAR  is  the  same  as  AVAR  for  white  and 
flicker  PM  or  white  FM  noise.  Bias  corrections  of  the  form  l/[l-a(x/T)],  where  T  is  the  record  length,  need  to  be 
applied  for  flicker  and  random  walk  FM  noise,  where  a  =  0.481  and  0.750,  respectively. 

The  number  of  equivalent  degrees  of  freedom  for  TOTVAR  can  be  estimated  for  white  FM,  flicker  FM  and 
random  walk  FM  noise  by  the  expression  b(T/T)— c,  where  b  =  1.500,  1.168  and  0.927,  and  c  =  0,  0.222  and  0.358, 
respectively.  For  white  and  flicker  PM  noise,  the  edf  for  a  total  deviation  estimate  is  the  same  as  that  for  the 
overlapping  ADEV  with  the  number  of     degrees  of  freedom  increased  by  2. 
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5.2.12.        Modified  Total  Variance 


The  modified  total  variance,  MTOT,  is  another  new  statistic  for  the  analysis  of 
frequency  stability.  It  is  similar  to  the  modified  Allan  variance,  MVAR,  and  has  the 
same  expected  value,  but  offers  improved  confidence  at  long  averaging  times.  It 
uses  the  same  phase  averaging  technique  as  MVAR  to  distinguish  between  white 
and  flicker  PM  noise  processes. 


The  modified  total  variance 
combines  the  features  of  the 
modified  Allan  and  total 
variances. 


A  calculation  of  MTOT  begins  with  an  array  of  N  phase  data  points  (time  deviates,  x,)  with  sampling  period  Tq  that 
are  to  be  analyzed  at  averaging  time  i  =  mto-  MTOT  is  computed  from  a  set  of  N  -  3m  +  1  subsequences  of  3m 
points.  First,  a  linear  trend  (frequency  offset)  is  removed  from  the  subsequence  by  averaging  the  first  and  last  halves 
of  the  subsequence  and  dividing  by  half  the  interval.  Then  the  offset-removed  subsequence  is  extended  at  both  ends 
by  uninverted,  even  reflection.  Next  the  modified  Allan  variance  is  computed  for  these  9m  points.  Finally,  these 
steps  are  repeated  for  each  of  the  N  -  3m  +  1  subsequences,  calculating  MTOT  as  their  overall  average.  These  steps, 
similar  to  those  for  MTOT.  but  acting  on  fractional  frequency  data,  are  shown  in  Figure  10. 


Phase  Data  x, ,  i  =  1  to  N 


N-3m+1  Subsequences; 


3m 


•  •  •       i=n  to  n+3m-1 


Linear  Trend  Removed: 

Extended  Subsequence- 
Unlnverted.  Even  Reflection 


9m 


X,  =  X,  -  c,  I,  c,  =  freq  offset 


0      »  0 

0  » 

0                »  _ 

•  •  •  ' 

1  f  ^ 

'   •  •  • 

Xn*3m-I  1  <  I  <  3m 


9  m-Point  Averages: 


6m  2nd  Differences: 


•  •  • 


Calculate  Mod  a,  (i)  for  Subsequence. 


mod  a,  ^(r)  =  1/2t  ■  (  z„  ^(m)  ),  where 
z,(m)  =  x„(m)  -  2x„,^(m)  +  x„,2^(m) 


Figure  10.   Steps  similar  to  calculation  of  MTOT  on  fractional  frequency  data. 


Computationally,  the  MTOT  process  requires  three  nested  loops: 

•  An  outer  summation  over  the  N  -  3m  +  1  subsequences.  The  3m-point  subsequence  is  formed,  its  linear  trend  is 
removed,  and  it  is  extended  at  both  ends  by  uninverted,  even  reflection  to  9m  points. 

•  An  inner  summation  over  the  6m  unique  groups  of  m-point  averages  from  which  all  possible  fully  overlapping 
second  differences  are  used  to  calculate  MVAR. 

•  A  loop  within  the  inner  summation  to  sum  the  phase  averages  for  three  sets  of  m  points. 

The  final  step  is  to  scale  the  result  according  to  the  sampling  period,  iq,  averaging  factor,  m,  and  number  of  points,  N. 
Overall,  this  can  be  expressed  as: 

1  N-im-l  f    1     n  +  lm-l 

ModTotvar(T)  =   Y     —  Y  {'z'(m)']'\,  (27) 

2(mrJ-(7V-3w  +  l)  tt   l6w,=t1„L  j  j 
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where  the  z,  (m)  terms  are  the  phase  averages  from  the  triply  extended  subsequence,  and  the  prefix  denotes  that  the 
linear  trend  has  been  removed.  At  the  largest  possible  averaging  factor,  m  =  N/3,  the  outer  summation  consists  of 
only  one  term,  but  the  inner  summation  has  6m  terms,  thus  providing  a  sizable  number  of  estimates  for  the  variance. 


Reference  for  Modified  Total  Variance 
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3 r'  PTTI  Meeting,  pp.  267-276,  Dec.  1999. 


5.2.13. 


Time  Total  Variance 


The  time  total  variance,  TTOT,  is  a  similar  measure  of  time  stability,  based  on  the 
modified  total  variance.  It  is  defined  as 


(28) 


The  time  total  variance  is  a 
measure  of  time  stability 
based  on  the  modified  total 
variance. 


5.2.14. 


Hadamard  Total  Variance 


The  Hadamard  total  variance,  HTOT,  is  a  total  version  of  the  Hadamard  variance. 
As  such,  it  rejects  linear  frequency  drift  while  offering  improved  confidence  at  large 
averaging  factors. 


An  HTOT  calculation  begins  with  an  array  of  N  fractional  frequency  data  points,  y, 
with  sampling  period  tq  that  are  to  be  analyzed  at  averaging  time  x  =m  Tq.  HTOT  is 
computed  from  a  set  of  N  -  3m  +  1  subsequences  of  3m  points.  First,  a  linear  trend 
(frequency  drift)  is  removed  from  the  subsequence  by  averaging  the  first  and  last 
halves  of  the  subsequence  and  dividing  by  half  the  interval.  Then  the  drift-removed 

subsequence  is  extended  at  both  ends  by  uninverted,  even  reflection.  Next  the  Hadamard  variance  is  computed  for 
these  9m  points.  Finally,  these  steps  are  repeated  for  each  of  the  N  -  3m  +  1  subsequences,  calculating  HTOT  as  their 
overall  average.  These  steps  are  shown  in  Figure!  1 . 


The  Hadamard  total  variance 
combines  the  features  of  the 
Hadamard  and  total  variances 
by  rejecting  linear  frequency 
drift,  handling  more  divergent 
noise  types,  and  providing 
better  confidence  at  large 
averaging  factors.  
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Fractional  Frequency  Data  y.  i  =  1  to  N 


N-3m+1  Subsequences 


•  •  • 


3m 


•      i=n  to  n+3m-1 


Linear  Freq  Drift  Removed: 

Extended  Subsequence 
Uninverted,  Even  Reflection 


9m 


y,  =  y,  -  q  i.  q  =  freq  drift 


0       »  0 

0  t 
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0  * 
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•  •  •  ' 

1  f  ' 

'    •  •  • 
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1  <  I <  3m 


9  m-Point  Averages: 


6m  2nd  Differences: 


Calculate  Hado   (t)  for  Subsequence 


Had  a,  ^t)  =  1/6  ■  (  z„  ^(m)  ),  where 


z„(m)  =  y„(m)  -  2y„,  Jm)  +  y„.2m(m) 
Then  Find  HTOT  as  Average  of  Subestimates 

Figure  11.  Steps  to  calculate  Hadamard  Total  Variance. 

Computationally,  the  HTOT  process  requires  three  nested  loops: 

•  An  outer  summation  over  the  N  -  3m  +  1  subsequences.  The  3m-point  subsequence  is  formed,  its  linear  trend  is 
removed,  and  it  is  extended  at  both  ends  by  uninverted,  even  reflection  to  9m  points. 

•  An  inner  summation  over  the  6m  unique  groups  of  m-point  averages  from  which  all  possible  fully  overlapping 
second  differences  are  used  to  calculate  HVAR. 

•  A  loop  within  the  inner  summation  to  sum  the  frequency  averages  for  three  sets  of  m  points. 

The  final  step  is  to  scale  the  result  according  to  the  sampling  period,  tq,  averaging  factor,  m,  and  number  of  points,  N. 
Overall,  this  can  be  expressed  as: 


1 


6{N-3m  +  \)  t 


/V-3m+\  (    1  n+3w-l 


(29) 


where  the  H|(m)  terms  are  the  Zn(m)  Hadamard  second  differences  from  the  triply  extended,  drift-removed 
subsequences.  At  the  largest  possible  averaging  factor,  m  =  N/3,  the  outer  summation  consists  of  only  one  term,  but 
the  inner  summation  has  6m  terms,  thus  providing  a  sizable  number  of  estimates  for  the  variance.  The  Hadamard 
total  variance  is  a  biased  estimator  of  the  Hadamard  variance,  so  a  bias  correction  is  required  that  is  dependent  on  the 
power  law  noise  type  and  number  of  samples. 

The  following  plots  shown  the  improvement  in  the  consistency  of  the  overlapping  Hadamard  deviation  results 
compared  with  the  normal  Hadamard  deviation,  and  the  extended  averaging  factor  range  provided  by  the  Hadamard 
total  deviation  [10]. 
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Figure  12.  (a)  Hadamard  Deviation,  (b)  Overlapping  Hadamard  Deviation. 
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Figure  13.  (a)  Hadamard  Total  Deviation,  (b)  Overlapping  &  Total  Hadamard  Deviations 


A  comparison  of  the  overlapping  and  total  Hadamard  deviations  shows  the  tighter  error  bars  of  the  latter,  allowing  an 
additional  point  to  be  shown  at  the  longest  averaging  factor. 


The  Hadamard  variance  may  also  be  used  to  perform  a  frequency  domain  (spectral)  analysis  because  it  has  a  transfer 
function  that  is  a  close  approximation  to  a  narrow  rectangle  of  spectral  width  1/(2  N  to),  where  N  is  the  number  of 
samples,  and  tq  is  the  measurement  time  [3].  This  leads  to  a  simple  expression  for  the  spectral  density  of  the 
fractional  frequency  fluctuations  Sy(f)  «  0.73  xo  Ho^yiT)  I  N,  where  f  =  1/  (2-io),  which  can  be  particularly  useful  at 
low  Fourier  frequencies. 

The  Picinbono  variance  is  a  similar  three-sample  statistic.  It  is  identical  to  the  Hadamard  variance  except  for  a  factor 
of  2/3  [4].  Sigma-z  is  another  statistic  that  is  similar  to  the  Hadamard  variance  that  has  been  applied  to  the  study  of 
pulsars  [5]. 
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It  is  necessary  to  identify  the  dominant  power  law  noise  type  as  the  first  step  in  determining  the  estimated  number  of 
chi-squared  degrees  of  freedom  for  the  Hadamard  statistics  so  their  confidence  limits  can  be  properly  set  [6].  Because 
the  Hadamard  variances  can  handle  the  divergent  flicker  walk  FM  and  random  run  FM  power  law  noises,  techniques 
for  those  noise  types  must  be  included.  Noise  identification  is  particularly  important  for  applying  the  bias  correction 
to  the  Hadamard  total  variance. 


References  for  Hadamard  Total  Variance 

1.  D.A.  Howe.  R.  Beard,  C.A.  Greenhall,  F.  Vernotte,  and  W.J.  Riley,  "A  total  estimator  of  the  Hadamard  function 
used  for  GPS  operations,"  Proc.  32"''  PTTI  Mtg.,  pp.  255-267,  Nov.  2000. 

2.  D.A.  Howe,  R.  Beard,  C.A.  Greenhall,  F.  Vernotte,  and  W.J.  Riley,  "Total  Hadamard  variance:  Application  to 
clock  steering  by  Kalman  filtering,"  Proc.  2001  European  Freq.  and  Time  Forum,  pp.  423-427  (Mar.  2001). 

3.  Chronos  Group,  Frequency  measurement  and  control.  Section  3.3.3,  Chapman  &  Hall,  London,  ISBN  0-412- 
48270-3  (1994). 

4.  B.  Picinbono,  «  Processus  a  accroissements  stationnaires  »  Ann.  des  telecom,  30(7-8):  211-212  (July- Aug. 
1975). 

5.  D.N.  Matsakis  and  F.J.  Jostles,  "Pulsar-appropriate  clock  statistics,"  Proc.  28'^  PTTI  Mtg.,  pp.  225-236  (Dec. 
1996). 

6.  D.A.  Howe,  R.L.  Beard,  C.A.  Greenhall,  F.  Vernotte,  W.J.  Riley,  and  T.K.  Peppier,  "Enhancements  to  GPS 
operations  and  clock  evaluations  using  a  total  Hadamard  deviation,"  IEEE  Trans.  Ultrason.  Ferroelect.  Freq. 
Cont.,  52:  1253-1261  (Aug.  2005). 


5.2.15.  Theo1 


Theol  is  a  new  class  of  statistics  which  mimic  the  properties  of  the  Allan  variance  Theol  is  a  two-sample  variance 

(AVAR)  while  covering  a  larger  range  of  averaging  times,  10  to  N-2  for  Theol  vs.  with  improved  confidence  and 

1  to  (N-l)/2  for  AVAR  [1].  It  provides  improved  confidence  and  the  ability  to  extended  averaging  factor 

obtain  a  result  for  a  maximum  averaging  time  equal  to  75  %  of  the  record  length.  range.  

Theol  [1]  is  defined  as  follows: 

I  N-m  m/  2-1  I  ^ 

where  m  =  averaging  factor,  tq-  measurement  interval,  and  jV=  number  of  phase  data  points,  for  m  even,  and  10  <  w 
<  N  -  \  .  It  consists  of  N  -  m  outer  sums  over  the  number  of  phase  data  points  -1,  and  m/2  inner  sums.  Theol  is  the 
rms  of  frequency  differences  averaged  over  an  averaging  time  i  =  0.75  (m  -  1)to. 

A  schematic  for  a  Theol  calculation  is  shown  in  Figure  14.  This  example  is  for  eleven  phase  samples  {N  =  1 1)  at  the 
largest  possible  averaging  factor  (m  =  10). 
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Theo1  Schematic  for  n=11 ,  m=  10 
i  =1  to  n-m  =  1  ,  5  =  0  to  m/2  -1=4 
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Figure  14.  A  schematic  for  Theol  calculation. 

The  single  outer  summation  (/  =  1  to  1)  at  the  largest  possible  averaging  factor  consists  of  ml2  =  5  ternis,  each  with 
two  phase  differences.  These  terms  are  scaled  by  their  spans  mil  ~  S  =  5  thru  1  so  that  they  all  have  equal  weighting. 
A  total  of  10  terms  contribute  to  the  Theol  statistic  at  this  largest-possible  averaging  factor.  The  averaging  time,  x, 
associated  with  a  Theol  value  is  t=  O.TS  m  io,  where  tq  is  the  measurement  interval.  Theol  has  the  same  expected 
value  as  the  Allan  variance  for  white  FM  noise,  but  provides  many  more  samples  that  provide  improved  confidence 
and  the  ability  to  obtain  a  result  for  a  maximum  t  equal  to  three-fourths  of  the  record  length,  T.  Theol  is  a  biased 
estimator  of  the  Allan  variance,  AVAR,  for  all  noise  types  except  white  FM  noise,  and  it  therefore  requires  the 
application  of  a  bias  correction.  Reference  [2]  contains  the  preferred  expression  for  determining  the  Theol  bias  as  a 
function  of  noise  type  and  averaging  factor: 

^1  ,  1  T^-       Avar  b 

Theol  Bias-  =  a+ — ,  .  (31) 

Theol  m'^ 

where  m  is  the  averaging  factor  and  the  constants  a,  b  and  c  are  given  in  Table  2.  Note  that  the  effective  tau  for  a 
Theol  estimation  is  t  =  O.TS  m  to,  where  to  is  the  measurement  interval. 


Table  2.  Theol  bias  parameters. 


Noise 

Alpha 

a 

b 

c 

RW  FM 

—2 

2.70 

-1.53 

0.85 

F  FM 

-1 

1.87 

-1.05 

0.79 

W  FM 

0 

1.00 

0.00 

0.00 

F  PM 

1 

0.14 

0.82 

0.30 

W  PM 

2 

0.09 

0.74 

0.40 

Empirical  formulae  have  been  developed  [1]  for  the  number  of  equivalent  x'  degrees  of  freedom  for  the  Theol 
statistic,  as  shown  in  Table  3: 
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Table  3.  Theol  EDF  Formulae 


Noise 

EDF 

RW  FM 

(  4.4  N-2  ]( (4.4    -  1)'  -  8.6r(4.4^  - 1)  + 1  \Ar) 

I    2.9r    )[                 (4.4iV-3)-  J 

F  FM 

r  2N'  -  l3Nr  -  3.5r  Y     r'  \ 

[           Nr  j[r'+2.3l 

W  FM 

"4.1A'+0.8    3.W  +  6.51 

r 

/ 

r               N  J 

,r^''^+5.2^ 

F  PM 

'4.798A^' -6.374A^r  +  12.387r  Y     r  \ 

(r  +  36.6f-(N-r)        jU  +  0.3j 

W  PM 

|^0.86(A^  +  l)(A^-4r/3)^ 

If     '  \ 

I  N-r 

where  r  =  0.75m,  and  with  the  condition  to<  T/10. 

5.2.16. 


TheoH 


Theol  has  the  same  expected  value  as  the  Allan  variance  if  bias  is  removed  [2].  It  is 
useful  to  combine  a  bias-removed  version  of  Theol,  called  TheoBR,  with  AVAR  to 
produce  a  composite  stabilit)'  plot.  The  composite  is  called  "TheoH"  which  is  short 
for  "hybrid-TheoBR"  [3].  TheoH  is  the  best  statistic  available  for  estimating  the 
stability  level  and  ty  pe  of  noise  of  a  frequency  source,  particularly  at  large  averaging 
times  and  with  a  mi.xture  of  noise  types  [4]. 


NewTheo  1 ,  TheoBR ,  and 
TheoH  are  versions  of  Theol 
that  provide  bias  removal  and 
combination  with  the  Allan 
variance. 


The  New  Theol  algorithm  of  Reference  [2]  provides  a  method  of  automatic  bias  correction  for  a  Theol  estimation 
based  on  the  average  ratio  of  the  Allan  and  Theol  variances  over  a  range  of  averaging  factors: 


NewTheo  1  (m,    ,  jV)= 


1 


■I- 


Avar(w=9+3i,ro,N) 


n+l'^Theol(m  =  12  +  4/.ro,^) 


Theol  (w,rQ,7V) 


(32) 


where  n  — 


N_ 
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,  and  |_  J  denotes  the  floor  function. 


NewTheo  1  was  used  in  Reference  [2]  to  form  a  composite  AVAR/  NewTheo  1  resuU  called  LONG,  which  has  been 
superseded  by  TheoH  (see  below). 

TheoBR  [3]  is  an  improved  bias-removed  version  of  Theol  given  by 


TheoBR(w,ro,iV) 


where  n  = 


3 


1   ^  Avar(m=9+3i,ro,N) 
n+ 1 Theo l(m  =  1 2  +  4/, To .  A^) 


,  and  |_  J  denotes  the  floor  function. 


Theol(m,ro,  A^), 


(33) 
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TheoBR  can  determine  an  unbiased  estimate  of  the  Allan  variance  over  the  widest  possible  range  of  averaging  times 
without  explicit  knowledge  of  the  noise  type.  TheoBR  in  the  equation  above  is  computationally  intensive  for  large 
data  sets,  but  computation  time  is  significantly  reduced  by  phase  averaging  with  negligible  effect  on  bias  removal  [5]. 

TheoH  is  a  hybrid  statistic  that  combines  TheoBR  and  AVAR  on  one  plot: 


Avar  (m.r,,  ^)  for  1  <  m  < 

TheoH(m,r„,iV)  .  ,, 

'  TheoBR  (m,r„>')  for 


0  75rn 


-<m<  N -\,m  even 


(34) 


where  k  is  the  largest  available  r  <  20%  T. 

An  example  of  a  TheoH  plot  is  shown  in  Figure  15: 
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Figure  15.  An  example  plot  of  TheoH. 


TheoH  is  a  composite  of  AVAR  and  bias-corrected  TheoBR  analysis  points  at  a  number  of  averaging  times 
sufficiently  large  to  form  a  quasi-continuous  curve.  The  data  are  a  set  of  1001  simulated  phase  values  measured  at 
15-minute  intervals  taken  over  a  period  of  about  10  days.  The  AVAR  results  are  able  to  characterize  the  stability  to 
an  averaging  time  of  about  two  days,  while  Theol  is  able  to  extend  the  analysis  out  to  nearly  a  week,  thus  providing 
significantly  more  information  from  the  same  data  set.  An  example  of  analysis  using  TheoH  with  data  from  a  Cs 
standard  is  shown  in  Section  1 1 . 
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References  for  Theo1,  NewTheol,  TheoBR  and  TheoH 

1.  D.A.  Howe  and  T.K.  Peppier,  "Very  Long-Term  Frequency  Stability:  Estimation  Using  a  Special-Purpose 
Statistic."  Proceedings  of  the  2003  IEEE  International  Frequency  Control  Symposium,  pp.  233-238,  May  2003. 

2.  D.A.  Howe  and  T.N.  Tasset,  "Theol:  Characterization  of  Very  Long-Term  Frequency  Stability,"  Proc.  2004 
EFTF. 

3.  D.A.  Howe,  "TheoH:  A  Hybrid,  High-Confidence  Statistic  that  Improves  on  the  Allan  Deviation,"  Metrologia  43 
(2006),  S322-S331. 

4.  J.  McGee  and  D.A.  Howe,  "TheoH  and  Allan  Deviation  as  Power-Law  Noise  Estimators,"  IEEE  Trans. 
Ultrasonics,  Ferroelectrics  and  Freq.  Contrl.,  Feb.  2007. 

5.  J.  McGee  and  D.A.  Howe,  "Fast  TheoBR:  A  method  for  long  data  set  stability  analysis,"  IEEE  Trans. 
Ultrasonics,  Ferroelectrics  and  Freq.  Contrl..  to  be  published,  2008. 


5.2.17  MTIE 


The  maximum  time  interval  error,  MTIE,  is  a  measure  of  the  maximum  time  error  of 
a  clock  over  a  particular  time  interval.  This  statistic  is  commonly  used  in  the 
telecommunications  industry.  It  is  calculated  by  moving  an  n-point  (n  =^  xIxq)  window 
through  the  phase  (time  error)  data  and  finding  the  difference  between  the  maximum 
and  minimum  values  (range)  at  each  window  position.  MTIE  is  the  overall 
maximum  of  this  time  interval  error  over  the  entire  data  set: 

MTIE{t)  =  A/ax„,,^_„  {^ax,,^<ic.n     )  -  (x.)},  (35) 

where  n  =  1,2,...,  N-1  and  N  =  number  of  phase  data  points. 

MTIE  is  a  measure  of  the  peak  time  deviation  of  a  clock  and  is  therefore  very  sensitive  to  a  single  extreme  value, 
transient  or  outlier.  The  time  required  for  an  MTIE  calculation  increases  geometrically  with  the  averaging  factor,  n, 
and  can  become  very  long  for  large  data  sets  (although  faster  algorithms  are  available  -  see  [1]  below). 

The  relationship  between  MTIE  and  Allan  variance  statistics  is  not  completely  defined,  but  has  been  the  subject  of 
recent  theoretical  work  [2,3].  Because  of  the  peak  nature  of  the  MTIE  statistic,  it  is  necessary  to  express  it  in  terms  of 
a  probability  level,  p,  that  a  certain  value  is  not  exceeded. 

For  the  case  of  white  FM  noise  (important  for  passive  atomic  clocks  such  as  the  most  common  rubidium  and  cesium 
frequency  standards),  MTIE  can  be  approximated  by  the  relationship 

MTIE(  r, /?)=k^  •  =  k^  •  V2  ■     ( r)  •  r ,  (36) 

where  kp  is  a  constant  determined  by  the  probability  level,  p,  as  given  in  Table  4,  and  ho  is  the  white  FM  power-law 
noise  coefficient. 


Table  4.  Constants  P  and  kp 


p,  % 

95 

1.77 

90 

1.59 

80 

1.39 

The  maximum  time  interval  error  (MTIE)  and  rms  time  interval  error  (TIE  rms)  are  clock  stability  measures 
commonly  used  in  the  telecom  industry  [4,  5].  MTIE  is  determined  by  the  extreme  time  deviations  within  a  sliding 


MTIE  IS  a  measure  of  clock 
error  commonly  used  in  the 
tele-communications 
industry.  
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window  of  span  x,  and  is  not  as  easily  related  to  such  clock  noise  processes  as  TDEV  [2].  MTIE  is  computationally 
intensive  for  large  data  sets  [6]. 

References  for  MTIE 

1.  S.  Bregni,  "Measurement  of  maximum  time  interval  error  for  telecommunications  clock  stability 
characterization,"  IEEE  Tram.  Instrum.  Meas.,  45(5):  900-906  (Oct.  1996). 

2.  P.  Travella  and  D.  Meo,  "The  range  covered  by  a  clock  error  in  the  case  of  white  FM,"  Proc.  30th  PTTI  Mtg., 
pp.  49-60  (Dec.  1998). 

3.  P.  Travella,  A.  Dodone,  and  S.  Leschiutta,  "The  range  covered  by  a  random  process  and  the  new  definition  of 
MTIE,"  Proc.  28th  PTTI  Mtg.,  pp.  1 19-123  (Dec.  1996). 

4.  S.  Bregni,  "Clock  Stability  Characterization  and  Measurement  in  Telecommunications,"  IEEE  Trans.  Instrum. 
Meas.A6{6):  1284-1294  (Dec.  1997). 

5.  G.  Zampetti,  "Synopsis  of  timing  measurement  techniques  used  in  telecommunications,"  Proc.  24th  PTTI  Mtg., 
pp.  313-326  (Dec.  1992). 

6.  S.  Bregni  and  S.  Maccabruni,  "Fast  computation  of  maximum  time  interval  error  by  binary  decomposition," 
IEEE  Trans.  Instrum.  Meas.,  49(6):  1240-1244  (Dec.  2000). 

7.  M.J.  Ivens,  "Simulating  the  Wander  accumulation  in  a  SDH  synchronisation  network,"  Master's  Thesis, 
University  College,  London,  U.K.  (^ov.  1997). 


5.2.18.  TIErms 


The  rms  time  interval  error,  TIE  rms.  is  another  clock  statistic  commonly  used  by  the  telecommunications  industry. 
TIE  rms  is  defined  by  the  expression 


where  n  =  1 ,2,...,  N-1  and  N  =  #  phase  data  points. 

For  no  frequency  offset,  TIE  rms  is  approximately  equal  to  the  standard  deviation  of  the  fractional  frequency 
fluctuations  multiplied  by  the  averaging  time.  It  is  therefore  similar  in  behavior  to  TDEV,  although  the  latter  properly 
identifies  divergent  noise  types. 


Reference  for  TIE  rms 

S.  Bregni,  "Clock  Stability  Characterization  and  Measurement  in  Telecommunications,"  IEEE  Trans.  Instrum. 
Meas.,  Vol.  46,  No.  6,  pp.  1284-1294,  Dec.  1997. 


5.2.19.        Integrated  Phase  Jitter  and  Residual  FM 

Integrated  phase  jitter  and  residual  FM  are  other  ways  of  expressing  the  net  phase  or  frequency  jitter  by  integrating  it 
over  a  certain  bandwidth.  These  can  be  calculated  from  the  amplitudes  of  the  various  power  law  terms. 

The  power  law  model  for  phase  noise  spectral  density  (see  section  6. 1 )  can  be  written  as 
SAf)=K-f\  (38) 
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where  S^'is  the  spectral  density  of  the  phase  fluctuations  in  rad"/Hz,/is  the  modulation  frequency,  K  is  amplitude  in 
rad",  and  x  is  the  power  law  exponent.  It  can  be  represented  as  a  straight  line  segment  on  a  plot  of  S(,{f)  in  dB  relative 
to  1  rad"^/Hz  versus  log/in  hertz.  Given  two  points  on  the  plot  (J],  dBi)  and /2,  dB:),  the  values  of  x  and  K  may  be 
determined  by 

dB,  -  dB, 

X  =  =          ,  (39) 

io-(iogy; -log/2) 

and 

^  =  10^'°  ^  (40) 

The  integrated  phase  jitter  can  then  be  found  over  this  frequency  interval  by 

Af-  =  f:S,{f)-df=t'K-f-df 

A^'=^{fr'-fr')foTx^-\  (41) 

x  +  l 

Af-  =K-(\ogf,-\ogf,)  forx  =  \. 

It  is  usually  expressed  as  A(j)  in  rms  radians. 

Similarly,  the  spectral  density  of  the  frequency  fluctuations  in  Hz"/Hz  is  given  by 

SXf)  =  K  ■  SAf)  =  f-  ■  Sjf)  =  K  ■        ,  (42) 

where  Vo  is  the  carrier  frequency  in  hertz,  and  Sjf)  is  the  spectral  density  of  the  fractional  frequency  fluctuations  (see 
Section  6. 1 ). 

The  integrated  frequency  jitter  or  residual  FM  is  therefore 


Af--t:sAf)-df=f^K-r'--df 


Af-=-^{fr'-fr)^orx^3  (43) 
A/-=^-(log/,-log/)forx  =  -3. 

It  is  usually  expressed  as  Af  in  rms  hertz. 

The  value  of  SJ^J)  in  dB  can  be  found  from  the  more  commonly  used  £(/)  measure  of  SSB  phase  noise  to  carrier  power 
ratio  in  dBc/Hz  by  adding  3  dB.  The  total  integrated  phase  noise  is  obtained  by  summing  the  Acj)'  contributions  from 
the  straight-line  approximations  for  each  power  law  noise  type.  The  ratio  of  total  phase  noise  to  signal  power  in  the 
given  integration  bandwidth  is  equal  to  10  log  A(j)~. 


1 
! 
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References  for  Integrated  Phase  Jitter  and  Residual  FM 

1.  W.  J.  Riley,  "Integrate  Phase  Noise  and  Obtain  Residual  FM,"  Microwaves,  August  1979,  pp.  78-80. 

2.  U.L.  Rohde,  Digital  PLL  Frequency  Synthesizers,  pp.  4 1 1  -4 1 8,  Prentice-Hall,  Englewood  Cliffs,  1 983. 

3.  D.A.  Howe  and  T.N.  Tasset,  "Clock  Jitter  Estimation  based  on  PM  Noise  Measurements,"  Proc.  2003  Joint  Mtg. 
IEEE  Intl.  Freq.  Cont.  Symp.  and  EFTF  Conf.,  pp.  541-546,  May  2003. 


5.2.20. 


Dynamic  Stability 


A  dynamic  stability  analysis  uses  a  sequence  of  sliding  time  windows  to  perform  a  dynamic  Allan  (DAVAR)  or 
Hadamard  (DHVAR)  analysis,  thereby  showing  changes  (nonstationarity)  in  clock  behavior  versus  time.  It  is  able  to 
detect  variations  in  clock  stability  (noise  bursts,  changes  in  noise  level  or  type,  etc.)  that  would  be  difficult  to  see  in  an 
ordinary  overall  stability  analysis.  The  results  of  a  dynamic  stability  analysis  are  presented  as  a  three-dimensional 
surface  plot  of  log  sigma  versus  log  tau  or  averaging  factor  as  a  function  of  time  or  window  number. 

An  example  of  a  DAVAR  plot  is  shown  below.  This  example  is  similar  to  the  one  of  Figure  2  in  Reference  [1], 
showing  a  source  with  white  PM  noise  that  changes  by  a  factor  of  2  at  the  middle  of  the  record. 
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Sigmo  Ronge: 
2,219e-02 
to 

2-495e  +  00 
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in  Clock  Noises  Using  the 
Dynamic  Allan  Variance" 


Analysis  Windows 
987  Windows 
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&  Step  Size  9 
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9.870e  +  02  sec 


DAVAR.DAT 


Figure  16.  Example  of  a  DAVAR  plot. 
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References  for  Dynamic  Stability 

1.  L.  Galleani  and  P.  Tavella,  "Tracking  Nonstationarities  in  Clock  Noises  Using  the  Dynamic  Allan  Variance," 
Proc.  2005  Joint  FCS/PTTI  Meeting. 

2.  L.  Galleani  and  P.  Tavella,  "The  Characterization  of  Clock  Behavior  with  the  Dynamic  Allan  Variance,"  Proc. 
2003  Joint  FCS/EFTF  Meeting,  pp.  239-244. 


5.3.  Confidence  Intervals 

It  is  wise  to  include  error  bars  (confidence  intervals)  on  a  stability  plot  to  indicate  the  degree  of  statistical  confidence 
in  the  numerical  results.  The  confidence  limits  of  a  variance  estimate  depend  on  the  variance  type,  the  number  of  data 
points  and  averaging  factor,  the  statistical  confidence  factor  desired,  and  the  type  of  noise.  This  section  describes  the 
use  of-//  statistics  for  setting  the  confidence  intervals  and  error  bars  of  a  stability  analysis. 

It  is  generally  insufficient  to  simply  calculate  a  stability  statistic  such  as  the  Allan  deviation,  thereby  finding  an 
estimate  of  its  expected  value.  That  determination  should  be  accompanied  by  an  indication  of  the  confidence  in  its 
value  as  expressed  by  the  upper  and  (possibly)  lower  limits  of  the  statistic  with  a  certain  confidence  factor.  For 
example,  if  the  estimated  value  of  the  Allan  deviation  is  1.0  x  10"",  depending  on  the  noise  type  and  size  of  the  data 
set,  we  could  state  with  95  %  confidence  that  the  actual  value  does  not  exceed  (say)  1.2  x  10"".  It  is  always  a  good 
idea  to  include  such  a  confidence  limit  in  reporting  a  statistical  result,  which  can  be  shown  as  an  upper  numeric  limit, 
upper  and  lower  numeric  bounds,  or  (equivalently)  error  bars  on  a  plot.  Even  though  those  confidence  limits  or  error 
bars  are  themselves  inexact,  they  should  be  included  to  indicate  the  validity  of  the  reported  result. 

If  you  are  unfamiliar  with  the  basics  of  confidence  limits,  it  is  recommended  that  an  introductory  statistics  book  be 
consulted  for  an  introduction  to  this  subject.  For  frequency  stability  analysis,  the  emphasis  is  on  various  variances, 
whose  confidence  limits  (variances  of  variances)  are  treated  with  chi-squared  (x')  statistics.  Strictly  speaking, 
statistics  apply  to  the  classical  standard  variance,  but  they  have  been  found  applicable  to  all  of  the  other  variances 
(Allan,  Hadamard,  total,  Theol,  etc.)  used  for  frequency  stability  analysis.  A  good  introduction  to  confidence  limits 
and  error  bars  for  the  Allan  variance  may  be  found  in  Reference  [1].  The  basic  idea  is  to  (1)  choose  an  single  or 
double-sided  confidence  limits  (upper  or  upper  and  lower  bounds),  (2)  choose  an  appropriate  confidence  factor  (e.g., 
95  %),  (3)  determine  the  number  of  equivalent  degrees  of  freedom  (edf),  (4)  use  the  inverse  yj'  distribution  to  find 
the  normalized  confidence  limit(s),  and  (5)  multiply  those  by  the  nominal  deviation  value  to  find  the  error  bar(s). 

5.3.1.         Simple  Confidence  Intervals 

The  simplest  confidence  interval  approximation,  with  no  consideration  of  the  noise  type,  sets  the  ±la  (68  %)  error 
bars  at  ±c7y(t)/VN,  where  N  is  the  number  of  frequency  data  points  used  to  calculate  the  Allan  deviation. 

A  more  accurate  determination  of  this  confidence  interval  can  be  made  by  considering  the  noise  type,  which  can  be 
estimated  by  the  Bl  bias  function  (the  ratio  of  the  standard  variance  to  the  Allan  variance).  That  noise  type  is  then  be 
used  to  determine  a  multiplicative  factor,  Kn,  to  apply  to  the  confidence  interval: 

Noise  Type  Kn 

Random  Walk  FM  0.75 

Flicker  FM  0.77 

White  FM  0.87 

Flicker  PM  0.99 

White  PM  0.99 
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5.3.2 


Chi-Squared  Confidence  Intervals 


Chi-squared  statistics  can  be  applied  to  calculate  single  and  double-sided  confidence  intervals  at  any  desired 
confidence  factor.  These  calculations  are  based  on  a  determination  of  the  number  of  degrees  of  freedom  for  the 
estimated  noise  type.  Most  stability  plots  show  ±la  error  bars  for  its  overlapping  Allan  deviation  plot. 

The  error  bars  for  the  modified  Allan  and  time  variances  are  also  determined  by  Chi-squared  statistics,  using  the 
number  of  MVAR  degrees  of  freedom  for  the  particular  noise  type,  averaging  factor,  and  number  of  data  points. 
During  the  Run  function,  noise  type  estimates  are  made  at  each  averaging  factor  (except  the  last,  where  the  noise  type 
of  the  previous  averaging  factor  is  used). 

Sample  variances  are  distributed  according  to  the  expression 

r=— (44) 

where  f-  is  the  Chi-square,  s^  is  the  sample  variance,  is  the  true  variance,  and  edf  is  the  equivalent  number  of 
degrees  of  freedom  (not  necessarily  an  integer).  The  edf  is  determined  by  the  number  of  analysis  points  and  the  noise 
type.  Procedures  exist  for  establishing  single-  or  double-sided  confidence  intervals  with  a  selectable  confidence 
factor,  based  on  statistics,  for  many  of  its  variance  functions.  The  general  procedure  is  to  choose  a  single-  or 
double-limited  confidence  factor,  p,  calculate  the  corresponding  value,  determine  the  edf  from  the  variance  type, 
noise  type  and  number  of  analysis  points,  and  thereby  set  the  statistical  limit(s)  on  the  variance.  For  double-sided 
limits, 

edf  ,    7        ,  edf 

and  cj-^^^  =  s-  -.  (45) 


rip, edf)  ^-(1 


5.4.  Degrees  of  Freedom 

The  equivalent  number  of  x'  degrees  of  freedom  (edf)  associated  with  a  statistical  variance  (or  deviation)  estimate 
depends  on  the  variance  type,  the  number  of  data  points,  and  the  type  of  noise  involved.  In  general,  the  progression 
from  the  original  two-sample  (Allan)  variance  to  the  overlapping,  total,  and  Theol  variances  has  provided  larger  edfs 
and  better  confidence.  The  noise  type  matters  because  it  determines  the  extent  that  the  points  are  correlated.  Highly 
correlated  data  have  a  smaller  edf  than  does  the  same  number  of  points  of  uncorrelated  (white)  noise.  An  edf 
determination  therefore  involves  (1)  choosing  the  appropriate  algorithm  for  the  particular  variance  type,  (2) 
determining  the  dominant  power  law  noise  type  of  the  data,  and  (3)  using  the  number  of  data  points  to  calculate  the 
corresponding  edf 

5.4.1 .         AVAR,  MVAR,  TVAR,  and  HVAR  EDF 

The  equivalent  number  of  x'  degrees  of  freedom  (edf)  for  the  Allan  variance  (AVAR),  the  modified  Allan  variance 
(MVAR)  and  the  related  time  variance  (TVAR),  and  the  Hadamard  variance  (HVAR)  is  found  by  a  combined 
algorithm  developed  by  C.A.  Greenhall,  based  on  its  generalized  autocovariance  function  [2]. 
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Figure  17.  Overlapping  ADEV  EDF  for  W  FM  Noise. 


100 


10 


^::::::::;::;:;:::::;::h:::::;i;;;:;h:::i;:::b::h:H::::::::;::::::::;i 

2  i  i  i 

1                    1  i 

0 

 ■  i  i  i--f --f  ^  [-  -1  t--i-i-f  ■ 

-2\XVi  i-  -Vi--  --|---f  --f  -i  --h 

 Full  EDF  Algorithm  vs  m  i- 

1  10  100 


Averaging  Factor,  m 
Figure  18.  Overlapping  ADEV  EDF  for  N=100. 
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Figure  19.  Overlapping  EDF  for  N-lOO  and  W  FM  Noise 

This  method  for  estimating  the  edf  for  the  Allan,  modified  Allan,  and  Hadamard  variances  supersedes  the  following 
somewhat  simpler  empirical  approximations  (which  may  still  be  used). 

The  equivalent  number  of  x~  degrees  of  freedom  (edf)  for  the  fully  overlapping  Allan  variance  (AVAR)  can  be 
estimated  by  the  following  approximation  formulae  for  each  power  law  noise  type: 


Table  5.  AVAR  approximation  formulae  for  each  power  law  noise  type. 


Power  law 
noise  type 

AVAR  edf,  where 
N  =  #  phase  data  points,  m  =  averaging  factor 

=  t/to 

W  PM 

(N  +  l){N-2m) 
2{N  -m) 

F  PM 

exp 

_   \    2m    )  \ 

'(2m  +  l)(7V-l)^T' 
4  JJ 

W  FM 

"3(iV-l)    2(A^-2)1  4nf 

2m           N  4w"+5 

F  FM 

2(^-2)' 
22>N  -  4.9 

For  m  =  1 

Am(N  +  'im) 

For  m  >1 

RW  FM 

Un-2)  1 

-3m(A^-l)  +  4m' 

m 
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The  edf  for  the  modified  Allan  variance  (MVAR)  can  be  estimated  by  the  same  expression  as  the  overlapping 
Hadamard  variance  (see  below)  with  the  arguments  changed  as  follows  (valid  for  -2  <  a  <  2):  MVAR  and  TVAR  edf 
for  N.  m  and  a  =  MVAR  edf  for  N  +  1 ,  m  and  a  -  2. 

The  edf  for  the  fulh  overlapping  Hadamard  variance  (HVAR)  can  be  found  by  an  earlier  algorithm  also  developed  by 
C.A.  Greenhall  based  on  its  generalized  autocovariance  function.  The  HVAR  edf  is  found  either  as  a  summation  (for 
small  m  cases  with  a  small  number  of  terms)  or  from  a  limiting  form  for  large  m,  where  1/edf  =  (l/p)(aO-al/p).  w  ith 
the  coefficients  as  follows: 


Table  6.  HVAR  edf  coefficients. 


Power  law 
noise  ty  pe 

HVAR  edf 
coefficients 

aO 

al 

W  FM 

7/9 

1/2 

F  FM 

1.00 

0.62 

RW  FM 

31/30 

17/28 

FW  FM 

1.06 

0.53 

RR  FM 

1.30 

0.54 

5.4.2.  TOTVAR  EDF 

The  edf  for  the  total  variance  (TOTVAR)  is  given  by  the  formula  b(T/T)  -  c,  where  T  is  the  length  of  the  data  record,  t 
is  the  averaging  time,  and  b  and  c  are  coefficients  that  depend  on  the  noise  type,  as  shown  in  Table  7; 


Table  7.  TOTVAR  edf  coefficients 


Power  law 
noise  tv  pe 

TOTVAR  edf 
Coefficients 

b 

c 

White  FM 

1.50 

0 

Flicker  FM 

1.17 

0.22 

Random  walk 
FM 

0.93 

0.36 

5.4.3.  MTOT  EDF 

The  edf  for  the  modified  total  variance  (MTOT)  is  given  by  the  same  formula  h{T/x)  -  c,  where  T  is  the  length  of  the 
data  record,  t  is  the  averaging  time,  and  b  and  c  are  coefficients  that  depend  on  the  noise  type  as  shown  in  Table  8: 


Table  8.  MTOT  edf  coefficients. 


Power  law 
noise  type 

MTOT  edf 
Coefficients 

b 

c 

White  PM 

1.90 

2.10 

Flicker  PM 

1.20 

1.40 

White  FM 

1.10 

1.20 

Flicker  FM 

0.85 

0.50 

Random  walk 
FM 

0.75 

0.31 

41 


5.4.4„         Theo1  /  TheoH  EDF 

The  equivalent  number  of  degrees  of  freedom  (edf)  for  the  Theol,  hence,  TheoBR  and  TheoH  variances,  is 
determined  by  the  following  approximation  formulae  for  each  power  low  noise  type. 


Table  9.  TheoBR  and  TheoH  approximation  formulae  for  each  power  law  noise  type. 


Power  law 

noi^ip  tvnp 

Theo  1  edf,  where 
N  =  #  phase  data  points,  t  =  0.75m,  m  =  averaging  factor  =  t/tq 

White  PM 

edf  = 

'0.86(A^,, +l)(A^,-f  r)^ 

(     ^  1 

-r 

V  J 

U  +  1-14  J 

Flicker  PM 

^4.798^; -6.374iV  r  +  12.387rV    r  ^ 

^df  =  TTT^    

■     ^        (r  +  36.6)"'(A^^-r)        )[t  +  0.3J 

White  FM 

edf  = 

'4.1A^^+0.8  3.l7V^+6.5' 

f  '''' 

\ 
/ 

T 

tr-+5.2 

Flicker  FM 

edf 

(2N]-\.?>N^j-3.5A(  ^ 

~[           Nj  JUY2.3J 

Random  walk 
FM 

( 4.4^  -2  Y (4.4^-1)'  -8.6r(4.47V  -l)  +  11.4rY 
edf  —   

■     I     2.9r     )[                 (4.4^,-3)-  J 
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5.5.   Noise  Identification 

Identification  of  the  dominant  power  law  noise  type  is  often  necessary  for  setting  confidence  intervals  and  making 
bias  corrections  during  a  frequency  stability  analysis.  The  most  effective  means  for  power  noise  identification  are 
based  on  the  Bi  and  R(n)  functions  and  the  lag  1  autocorrelation. 

5.5.1.  Power  Law  Noise  Identification 

It  is  often  necessary  to  identify  the  dominant  power  law  noise  process  (WPM,  FPM,  WFM,  FFM,  RWFM,  FWFM  or 
RRFM)  of  the  spectral  density  of  the  fractional  frequency  fluctuations,  Sy{f)  =  hj  °-  (a  =  2  to  -4),  to  perform  a 
frequency  stability  analysis.  For  example,  knowledge  of  the  noise  type  is  necessary  to  determine  the  equivalent 
number  of  chi-squared  degrees  of  freedom  (edf)  for  setting  confidence  intervals  and  error  bars,  and  it  is  essential  to 
know  the  dominant  noise  type  to  correct  for  bias  in  the  newer  Total  and  Theol  variances.  While  the  noise  type  may 
be  known  a  priori  or  estimated  manually,  it  is  desirable  to  have  an  analytic  method  for  power  law  noise  identification 
that  can  be  used  automatically  as  part  of  a  stability  analysis  algorithm. 
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There  is  little  literature  on  the  subject  of  power-law  noise  identification.  The  most  common  method  for  power  law 
noise  identification  is  simply  to  observe  the  slope  of  a  log-log  plot  of  the  Allan  or  modified  Allan  deviation  versus 
averaging  time,  either  manually  or  by  fitting  a  line  to  it.  This  obviously  requires  at  least  two  stability  points.  During 
a  stability  calculation,  it  is  desirable  (or  necessary)  to  automatically  identify  the  power  law  noise  type  at  each  point, 
particularly  if  bias  corrections  and/or  error  bars  must  be  applied. 

5.5.2.  Noise  Identification  Using  Bi  and  R(n) 

A  noise  identification  algorithm  that  has  been  found  effective  in  actual  practice,  and  that  works  for  a  single  t  point 
over  the  full  range  of  -4  <  a  <  2  is  based  on  the  Barnes  Bi  function,  which  is  the  ratio  of  the  N-sample  (standard) 
variance  to  the  two-sample  (Allan)  variance,  and  the  R(n)  function  [1],  which  is  the  ratio  of  the  modified  Allan  to  the 
normal  Allan  variances.  The  B|  function  has  as  arguments  the  number  of  frequency  data  points,  N,  the  dead  time 
ratio,  r  (which  is  set  to  1),  and  the  power  law  t-domain  exponent,  |i.  The  Bi  dependence  on  (i  is  used  to  determine  the 
power  law  noise  type  for  -2  <  n  <  2  (W  and  F  PM  to  FW  FM).  For  a  Bi  corresponding  to  ^  =  -2,  the  a  =  1  or  2  (F  PM 
or  W  PM  noise)  ambiguity  can  be  resolved  with  the  R(n)  ratio  using  the  modified  Allan  variance.  For  the  Hadamard 
variance,  for  which  RR  FM  noise  can  apply,  (m  =3,  a  =  -4),  the  Bi  ratio  can  be  applied  to  frequency  (rather  than 
phase)  data,  and  adding  2  to  the  resulting  ^. 

The  overall  noise  B|/R(n)  noise  identification  process  is  therefore: 

1 .  Calculate  the  standard  and  Allan  variances  for  the  applicable  i  averaging  factor. 

/V(l-/V") 

2.  Calculate  B|(N,  r=l,  n)=  . 

2(/V  -  l)(l  -2") 

3.  Determine  the  expected  B|  ratios  for  a  =  -3  through  1  or  2. 

4.  Set  boundaries  between  them  and  find  the  best  power  law  noise  match. 

5.  Resolve  an  a  =  1  or  2  ambiguity  with  the  modified  Allan  variance  and  R(n). 

6.  Resolve  an  a  =  -3  or  ^  ambiguity  by  applying  B]  to  frequency  data. 

The  boundaries  between  the  noise  types  are  generally  set  as  the  geometric  means  of  their  expected  values.  This 
method  cannot  distinguish  between  W  and  F  PM  at  unity  averaging  factor. 

5.5.3.  The  Autocorrelation  Function 

The  autocorrelation  function  (ACF)  is  a  fundamental  way  to  describe  a  time  series  by  multiplying  it  by  a  delayed 
version  of  itself,  thereby  showing  the  degree  by  which  its  value  at  one  time  is  similar  to  its  value  at  a  certain  later 
time.  More  specifically,  the  autocorrelation  at  lag  k  is  defined  as 

_  E[{z,-^)iz,^,-^i)] 
Pk  1   '  (46) 

where  z,  is  the  time  series,  \jl  is  its  mean  value,  is  its  variance,  and  E  denotes  the  expected  value.  The 
autocorrelation  is  usually  estimated  by  the  expression 


2  N-k 


k  ■,     N  ' 

where  z  is  the  mean  value  of  the  time  series  and     is  the  number  of  data  points  [2]. 


(47) 
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5.5.4.  The  Lag  1  Autocorrelation 


The  lag  1  autocorrelation  is  simply  the  value  of  r\  as  given  by  the  expression  above.  For  frequency  data,  the  lag  1 
autocorrelation  is  able  to  easily  identify  white  and  flicker  PM  noise,  and  white  (uncorrelated)  FM  noise,  for  which  the 
expected  values  are  -1/2,  -1/3  and  zero,  respectively.  The  more  divergent  noises  have  positive  r\  values  that  depend 
on  the  number  of  samples,  and  tend  to  be  larger  (approaching  1).  For  those  more  divergent  noises,  the  data  are 
differenced  until  they  become  stationary,  and  the  same  criteria  as  for  WPM,  FPM  and  WFM  are  then  used,  corrected 
for  the  differencing.  The  results  can  be  rounded  to  determine  the  dominant  noise  type  or  used  directly  to  estimate  the 
noise  mixture. 

5.5.5.  Noise  Identification  Using 

An  effective  method  for  identifying  power  law  noises  using  the  lag  1  autocorrelation  [3]  is  based  on  the  properties  of 
discrete-time  fractionally  integrated  noises  having  spectral  densities  of  the  form  (2  sin  k For  5  <  Vi,  the  process 
is  stationary  and  has  a  lag  1  autocorrelation  equal  to  pi  =  5  /  (1-5)  [4],  and  the  noise  type  can  therefore  be  estimated 
from  5  =  ri  /  (1+ri).  For  frequency  data,  white  PM  noise  has  pi  =  -1/2,  flicker  PM  noise  has  pi  =  -1/3,  and  white  FM 
noise  has  pi  =  0.  For  the  more  divergent  noises,  first  differences  of  the  data  are  taken  until  a  stationary  process  is 
obtained  as  determined  by  the  criterion  5  <  0.25.  The  noise  identification  method  therefore  uses  p  =  -round  (25)  -2d, 
w  here  round  (25)  is  25  rounded  to  the  nearest  integer  and  d  is  the  number  of  times  that  the  data  is  differenced  to  bring 
5  down  to  <  0.25.  If  z  is  a  x-average  of  frequency  data  y(t),  then  a  =  p;  if  z  is  a  i-sample  of  phase  data  x(t),  then  a  =  p 
+  2,  where  a  is  the  usual  power  law  exponent  /°,  thereby  determining  the  noise  type  at  that  averaging  time.  The 
properties  of  this  power  law  noise  identification  method  are  summarized  in  Table  10.  It  has  excellent  discrimination 
for  all  common  power  law  noises  for  both  phase  and  frequency  data,  including  difficult  cases  with  mixed  noises. 


Noise 
Type 

a 

Phase  Data* 
x(t) 

d=0  ACF  of 
Phase  Data 

W  PM 

2 

F  PM 

1 

m 

V 

W  FM 

0 

/ 

F  FM 

-1 

A 

RW  FM 

-2 

*  The  differencing  operation  changes  the  appearance 
of  the  phase  data  to  that  shown  2  rows  higher. 

Figure  20.  Lag  1  Autocorrelation  for  Various  Power  Law  Noises 
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Table  10.  Lag  1  Autocorrelation  for  Various  Power  Law  Noises  and  Differences 


Noise 
Type 

a 

Lag  1  Autocorrelation,  r/ 

d=0 

d=l 

d=2 

x(t) 

y(t) 

x(t) 

y(t) 

x(t) 

y(t) 

~} 

0 

-1/2 

-1/2 

-2/3 

-2/3 

-3/4 

1 

=0.7 

-1/3 

-1/3 

-3/5 

-3/5 

-5/7 

0 

=  1 

0 

0 

-1/2 

-1/2 

-2/3 

-1 

=  1 

=0.7 

=0.7 

-1/3 

-1/3 

-3/5 

.2 

*1 

=  1 

=  1 

0 

0 

-1/2 

Shaded  values  are  those  used  for  noise  ID  for  the  particular 
noise  and  data  type. 

5.5.6.  Noise  ID  Algorithm 

The  basic  lag  1  autocorrelation  power  law  noise  identification  algorithm  is  quite  simple.  The  inputs  are  a  vector  zj,..., 
-N  of  phase  or  fi-equency  data,  the  minimum  order  of  differencing  dmin  (default  =  0),  and  the  maximum  order  of 
differencing  dmax.  The  output  is  p,  an  estimate  of  the  a  of  the  dominant  power  law  noise  type,  and  (optionally)  the 
value  of  d. 


Done  =  False,  d  =  0 
While  Not  Done 

1  ^ 
z  = — Vz, 



1  +  r, 

If  d  >=  dmin  And  (5  <  0.25  Or  d  >=  dmax) 
p  =  -2iS  +  d) 
Done  =  True 

Else 

^1  —  ^2  ~  ■^i'---'-^A'-i  ~      ~  ■^yv-i 
TV  =     - 1 

d  =  d  +  \ 

End  If 
End  While 

Note:  May  round p  to  nearest  integer 

Figure  21.  The  basic  lag  1  autocorrelation  power  law  noise  identification  algorithm. 
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The  input  data  should  be  for  the  particular  averaging  time,  t,  of  interest,  and  it  may  therefore  be  necessary  to  decimate 
the  phase  data  or  average  the  frequency  by  the  appropriate  averaging  factor  before  applying  the  noise  identification 
algorithm.  The  dmax  parameter  should  be  set  to  two  or  three  for  an  Allan  or  Hadamard  (two  or  three-sample) 
variance  analysis,  respectively.  The  alpha  result  is  equal  to  p+2  or  p  for  phase  or  frequency  data,  respectively,  and 
may  be  rounded  to  an  integer  (although  the  fractional  part  is  useful  for  estimated  mixed  noises).  The  algorithm  is  fast, 
requiring  only  the  calculation  of  one  autocorrelation  value  and  first  differences  for  several  times.  It  is  independent  of 
an\  particular  variance.  The  lag  1  autocorrelation  method  yields  good  results,  consistently  identifying  pure  power 
noise  for  a  =  2  to  ^  for  sample  sizes  of  about  30  or  more,  and  generally  identifying  the  dominant  type  of  mixed 
noises  when  it  is  at  least  10  %  larger  than  the  others.  For  a  mixture  of  adjacent  noises,  the  fractional  result  provides 
an  indication  of  their  ratio.  It  can  handle  all  averaging  factors. 

Before  analysis,  the  data  should  be  preprocessed  to  remove  outliers,  discontinuities,  and  deterministic  components. 
Acceptable  results  can  be  obtained  from  the  lag  1  autocorrelation  noise  identification  method  for  N  >  32,  where  N  is 
the  number  of  data  points.  The  algorithm  tends  to  produce  jumps  in  the  estimated  alpha  for  mixed  noises  when  the 
differencing  factor,  d,  changes  (although  the  alpha  value  when  rounded  to  an  integer  is  still  consistent).  This  can  be 
avoided  by  using  the  same  d  for  the  entire  range  of  averaging  times,  at  the  expense  of  higher  variability  when  a  lower 
d  would  have  been  sufficient.  The  lag  1  autocorrelation  method  for  power  law  noise  identification  is  a  fast  and 
effective  way  to  support  the  setting  of  confidence  intervals  and  to  apply  bias  corrections  during  a  frequency  stability 
analysis,  as  shown  Figure  22: 


SAO  VLG1 1 B  H-Maser  S/N  SA026  vs  SA01 8 

Noise  Type 


Averaging  Time.t,  Seconds 


Figure  22.  Frequency  Stability  and  Noise  Analysis  of  Two  Hydrogen  Masers 
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5.6.   Bias  Functions 

Several  bias  functions  are  defined  and  used  in  the  analysis  of  frequency  stability,  as  defined  below.  In  particular,  B], 
the  ratio  of  the  standard  variance  to  the  Allan  variance,  and  R(n),  the  ratio  of  the  modified  Allan  variance  to  the 
normal  Allan  variance,  are  used  for  the  identification  of  power  law  noise  types  (see  section  5.2.2),  and  the  B2  and  B3 
bias  functions  are  used  to  correct  for  dead  time  in  a  frequency  stability  measurement. 


The  Bl  bias  function  is  the  ratio  of  the  N-sample  (standard)  variance  to  the  two-sample  (Allan)  variance  with  dead 
time  ratio  r  =  T/t,  where  T  =  time  between  measurements,  i  =  averaging  time,  and  \i  =  exponent  of  t  in  Allan  variance 
for  a  certain  power  law  noise  process: 


The  Bl  bias  fiinction  is  useful  for  performing  power  law  noise  identification  by  comparing  its  actual  value  to  those 
expected  for  the  various  noise  types  (see  section  5.2.2). 

5.8.   B2  Bias  Function 

The  Bt  bias  function  is  the  ratio  of  the  two-sample  (Allan)  variance  with  dead  time  ratio  r  =  T/t  to  the  two-sample 
(Allan)  variance  without  dead  time  (r  =  1 ): 


5.7.   Bi  Bias  Function 


B,(7V,r,//) 


cj'{N,T,T) 
(j\2,T,T)  ■ 


(48) 


(t\2,T,t) 
cr\2,T,T) 


(49) 
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5.9.   Bs  Bias  Function 

The  B3  bias  function  is  the  ratio  of  the  N-sample  (standard)  variance  with  dead  time  ratio  r  =  T/t  at  multiples  M  =  i/t, 
of  the  basic  averaging  time  tq  to  the  N-sample  variance  with  the  same  dead  time  ratio  at  averaging  time  x: 

B,{N,M.r.^)^     \/   1  \  '    .  (50) 
a-(N,T,T) 


The  product  of  the  Bi  and  B3  bias  functions  is  used  for  dead  time  correction,  as  discussed  in  section  5.7. 

5.10.  R(n)  Bias  Function 

The  R(n)  function  is  the  ratio  of  the  modified  Allan  variance  to  the  normal  Allan  variance  for  n  =  number  of  phase 
data  points.  Note:  R(n)  is  also  a  function  of  a,  the  exponent  of  the  power  law  noise  type: 

Mod  cr"(r) 

R{n)  =  .  (51) 

The  R(n)  bias  fiinction  is  useful  for  performing  power  law  noise  identification  by  comparing  its  actual  value  to  those 
expected  for  the  various  noise  ty  pes  (see  Section  5.2.2). 


5.11.  TOTVAR  Bias  Function 


The  TOTVAR  statistic  is  an  unbiased  estimator  of  the  Allan  variance  for  white  and  flicker  PM  noise,  and  for  white 
FM  noise.  For  flicker  and  random  walk  FM  noise,  TOTVAR  is  biased  low  as  i  becomes  significant  compared  with 
the  record  length.  The  ratio  of  the  expected  value  of  TOTVAR  to  AVAR  is  given  by  the  expression 


B(TOTAL)=l-a 


T 

0<r<  — ,  (52) 
2 


where  a  =  l/31n2  =  0.481  for  flicker  FM  noise,  a  =  3/4  =  0.750  for  random  walk  FM  noise,  and  T  is  the  record  length. 
At  the  maximum  allowable  value  of  x  =  T/2,  TOTVAR  is  biased  low  by  about  24  %  for  RW  FM  noise.  This  bias 
function  should  be  used  to  correct  all  reported  TOTVAR  results. 


5.12.  MTOT  Bias  Function 


The  MTOT  statistic  is  a  biased  estimator  of  the  modified  Allan  variance.  The  MTOT  bias  factor  (the  ratio  of  the 
expected  value  of  Mod  Totvar  to  MVAR),  as  shown  in  Table  11,  depends  on  the  noise  type  but  is  essentially 
independent  of  the  averaging  factor  and  number  of  data  points. 


Table  11.  MTOT  bias  factors  for  each  noise  type. 


Noise 

Bias  Factor 

W  PM 

1.06 

F  PM 

1.17 

W  FM 

1.27 

F  FM 

1.30 

RW  FM 

1.31 
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This  bias  factor  should  be  used  to  correct  all  reported  MTOT  results. 

5.13.  Theo1  Bias 

The  Theol  statistic  is  a  biased  estimator  of  the  Allan  variance.  The  Theol  bias  factor  (the  ratio  of  the  expected  value 
of  Theol  to  AVAR)  depends  on  both  noise  type  and  averaging  factor: 

^,  .  ,  r.-        AVAR  b 

Theol  Bias  =  =  a  +  — -,  (53) 

Theol 

where  m  is  the  averaging  factor  and  the  constants  a,  b  and  c  are  given  in  Table  12.  Note  that  the  effective  tau  for  a 
Theol  estimation  is  t  =  0.75-m-to,  where  to  is  the  measurement  interval. 


Table  12.  Constants  a,  b,  and  c  for  Theol  bias. 


Noise 

Alpha 

a 

b 

c 

RW  FM 

_2 

2.70 

-1.53 

0.85 

F  FM 

-1 

1.87 

-1.05 

0.79 

W  FM 

0 

1.00 

0.00 

0.00 

F  PM 

1 

0.14 

0.82 

0.30 

W  PM 

2 

0.09 

0.74 

0.40 

5.14.  TheoH  Bias 

TheoH  statistic  is  a  bias-removed  estimator  of  the  Allan  variance.  TheoH  is  the  best  statistic  for  estimating  the 
frequency  stability  and  type  of  noise  over  the  widest  possible  range  of  averaging  times  without  explicit  knowledge  of 
the  noise  type. 
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5.15.  Dead  Time 


Dead  time  can  occur  in  frequency  measurements  because  of 
instrumentation  delay  between  successive  measurements,  or  because  of  a 
deliberate  wait  between  measurements.  It  can  have  a  significant  effect  on 
the  results  of  a  stability  analysis,  especially  for  the  case  of  large  dead  time 
(e.g.,  frequency  data  taken  for  100  seconds,  once  per  hour). 


Dead  time  can  occur  in  frequency 
measurements  and  can  significantly 
affect  a  subsequent  stability  analysis. 
Methods  are  available  to  correct  for  dead 
time  and  thus  obtain  unbiased  results. 


Measurement 
Time 


Dead 
Time 


Measurement 
Time 


Dead 
Time 


Figure  23.  Illustration  of  dead  time  between  successive  measurements. 

Dead  time  corrections  can  be  applied  by  dividing  the  calculated  Allan  deviation  by  the  square  root  of  the  product  of 
the  Barnes  B^  and  B3  bias  ratios.  These  corrections  are  particularly  important  for  non-white  FM  noise  with  a  large 
dead  time  ratio.  Restricting  the  dead  time  corrections  to  Allan  deviations  is  a  conservative  approach  based  on  the  B2 
and  B3  definitions.  Those  bias  functions  depend  critically  on  the  power  law  noise  type.  Requiring  manual  noise 
selection  avoids  the  problem  of  noise  identification  for  biased  data  having  the  wrong  sigma-tau  slope.   Dead  time 
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correction  is  problematic  for  data  having  multiple  noise  types.  In  addition  to  introducing  bias,  measurement  dead 
time  reduces  the  confidence  in  the  results,  lowers  the  maximum  allowable  averaging  factor,  and  prevents  proper 
conversion  of  frequency  to  phase.  Moreover,  no  information  is  available  about  the  behavior  of  the  device  under  test 
during  the  dead  time.  It  is  recommended  that  these  issues  be  avoided  by  making  measurements  with  zero  dead  time. 

Dead  time  that  occurs  at  the  end  of  a  measurement  can  be  adjusted  for  in  an  Allan  deviation  determination  by  using 
the  Barnes  Bi  bias  function  [1],  the  ratio  of  the  two-sample  variance  with  dead  time  ratio  r  =  T/x  to  the  two-sample 
variance  without  dead  time.  Otherwise,  without  this  correction,  one  can  determine  only  the  average  frequency  and  its 
drift.  When  such  data  are  used  to  form  frequency  averages  at  longer  tau,  it  is  necessary  to  also  use  the  B3  bias 
function  [2],  the  ratio  of  the  variance  with  distributed  dead  time  to  the  variance  with  all  the  dead  time  at  the  end. 
Those  bias  corrections  are  made  by  use  of  the  product  of  Bt  and  B3.  The  power  law  noise  type  must  be  known  in 
order  to  calculate  these  bias  functions.  Simulated  periodically  sampled  frequency  data  with  distributed  dead  time  for 
various  power  law  noise  processes  shows  good  agreement  with  the  B2  and  B3  bias  function  corrections,  as  shown  in 
Figure  24. 
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Figure  24.  Frequency  stability  plots  for  common  power  law  noises  with  large  measurement  dead  time  (r  =  T/t  =  36). 
Simulated  data  sampled  for  x  =  100  seconds  once  per  hour  for  10  days.  Nominal  1x10  "  stability  at  x  =  100  seconds  shown 
by  lines.  Plots  show  stability  of  simulated  data  sets  for  continuous,  sampled  and  dead  time-corrected  data,  (a)  White  PM 
(^i  =  -2)  Vb.  =  0.82  at  AF  =  1,  (b)  Flicker  PM     =  -2),  VBj  =  0.82  at  AF  =  1,  (c)  White  FM     -  -1),  Vb,  =  1.00  at  AF  =  1,  (d) 
Flicker  FM  (|i  =  0)  Vb.  =  1.92  at  AF  =  1,  (e)  RW  FM  (|i  =  1)  VBj  =  7.31  at  AF  =  1. 
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These  simulations  show  that  the  B2  and  B3  bias  corrections  are  able  to  support  reasonably  accurate  results  for  sampled 
frequency  stability  data  having  a  large  dead  time,  when  the  power  law  noise  type  is  known.  The  slope  of  the  raw 
sampled  stability  plot  does  not  readily  identify  the  noise  type,  however,  and  mixed  noise  types  would  make  correction 
difficult.  The  relatively  small  number  of  data  points  reduces  the  confidence  of  the  results,  and  limits  the  allowable 
averaging  factor.  Moreover,  the  infrequent  measurements  provide  no  information  about  the  behavior  of  the  clock 
during  the  dead  time,  and  prevent  a  proper  conversion  of  frequency  to  phase.  Sparsely  sampled  data  are  therefore  not 
recommended  for  the  purpose  of  stability  analysis. 
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5.16.  Unevenly  Spaced  Data 

Unevenly  spaced  phase  data  can  be  handled  if  they  have  associated  timetags  by  using  the  individual  timetag  spacing 
when  converting  them  to  frequency  data.  Then,  if  the  tau  differences  are  reasonably  small,  the  data  may  be  analyzed 
by  use  of  the  average  timetag  spacing  as  the  analysis  tau,  in  effect  placing  the  frequency  data  on  an  average  uniform 
grid.  While  completely  random  data  spacing  is  not  amenable  to  this  process,  tau  variations  of  ±10  %  will  yield 
reasonable  results  as  long  as  the  exact  interval  is  used  for  phase  to  frequency  conversion. 

An  example  of  unevenly  spaced  data  is  two-way  satellite  time  and  frequency  transfer  (TWSTFT)  measurements  made 
on  Monday,  Wednesday,  and  Friday  of  each  week,  where  the  data  spacing  is  either  one  or  two  days. 


t>a(«:  09/29/03    Titr^:  09:i9-.i5 


FREQUENCY  STABILITY 


Composite  TDEV  Plot 


10'      1.5  2 
Averaging  Time,  T,  Days 


Figure  25.  TDEV  results  for  simulated  TWSTFT  data. 

The  TWSTFT  data  are  simulated  as  256  points  of  white  PM  noise  with  an  Allan  deviation  (ADEV)  level  of  ay(T)  =  I 
xlO  "  at  1-day.  A  composite  plot  of  the  TWSTFT  TDEV  results  is  shown  above.  The  corresponding  TDEV  is  5.77 


xlO"'-  sec  at  T  =  1  day  (TDEV 


MDEV  divided  by  V3),  as  shown  in  curve  A.  Note  that  these  time-stability  plots 
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include  points  at  all  possible  tau  values.  The  green  line  shows  that  the  -0.5  slope  of  the  TDEV  plot  for  W  PM  noise. 
The  TWSTFT  data  are  sampled  once  on  Monday.  Wednesday,  and  Friday  of  each  week.  These  sampled  data 
therefore  have  an  average  tau  of  7/3  =  2.33  days,  and  their  TDEV  is  shown  in  curve  B.  If  the  missing  points  are 
replaced  by  linearh  interpolated  phase  values,  the  TDEV  becomes  highly  distorted,  as  shown  in  curve  C.  If  the 
sampled  phase  data  are  converted  to  frequency  data  using  their  timetag  differences  to  determine  the  individual 
measurement  intervals,  the  average  tau.  i^.g.  is  close  to  2.33  days  (depending  on  the  final  fractional  week  that  is 
included),  and  the  resulting  TDEV  is  shown  in  curve  D.  It  is  essentially  identical  to  that  for  the  sampled  phase  data 
shown  in  curve  B.  It  is  interesting  that,  although  the  converted  frequency  results  are  different  depending  on  whether 
the  average  or  individual  (more  correct)  taus  are  used,  the  (integrated)  TDEV  results  are  not  (at  least  for  a  white  PM 
noise  process). 

None  of  the  results  is  in  good  agreement  with  the  nominal  simulation.  The  result  with  the  linearly  interpolated  phase 
points  is  particular!)  bad  for  Kia^g.  and  is  similar  to  that  of  Tavella  and  Leonardi.  as  shown  in  Figure  1  of  Reference 
[1].  As  they  point  out  in  that  paper,  because  the  true  sampling  interval  is  lave.  it  is  not  possible  to  estimate  the  noise  at 
shorter  times,  especialh  for  an  uncorrelated  white  noise  process.  They  further  suggest  that  the  higher  level  of  the 
estimated  noise  is  related  to  the  ratio  of  the  true  and  interpolated  sampling  times  (~2.33)  and  the  Vx  dependence  of 
TDEV.  By  applying  a  correction  factor  of  V2.33  1.5,  the  longer-tau  TDEV  estimates  are  lowered  to  the  correct 
level.  These  factors  are  smaller  for  other  non-white  PM  and  FM  noise  processes.  The  adjusted  method  of  using 
frequency  data  converted  from  phase  data  by  using  individual  tau  values  adjusted  for  the  timetag  spacing  is 
recommended  because  it  does  not  use  interpolation,  does  not  present  results  at  unrealisticall>  low  tau.  and  uses  the 
best  frequency  estimates. 

Another  situation  is  data  that  are  taken  in  bursts.  In  that  case,  the  best  approach  is  probably  to  analyze  the  segments 
separately,  perhaps  averaging  those  results  to  obtain  better  statistical  confidence.  One  could  obtain  reasonable  results 
for  the  shorter  av  eraging  times,  but  cannot  apply  standard  techniques  to  analyze  the  complete  data  set. 


References  for  Unevenly  Spaced  Data 

1.  P.  Tavella  and  M.  Leonardi.  "Noise  Characterization  of  Irregularly  Spaced  Data,"  Proceedings  of  the  12'^ 
European  Frequency  and  Time  Forum,  pp.  209-214.  March  1998. 

2.  C.  Hackman  and  T.E.  Parker.  "Noise  Analv  sis  of  Unevenly  Spaced  Time  Series  Data."  Metrologia.  Vol.  33.  pp. 
457-466.  1996. 


5.17.  Histograms 

A  histogram  shows  the  amplitude  distribution  of  the  phase  or  frequency  fluctuations,  and  can  provide  insight 
regarding  them.  We  can  expect  a  normal  (Gaussian)  distribution  for  a  reasonably  sized  data  set,  and  a  different  (e.g., 
bimodal)  distribution  can  be  a  sign  of  a  problem. 

For  a  normal  distribution,  the  standard  deviation  is  approximately  equal  to  the  half-width  at  half-height  (HWHA  = 
1.177s) . 
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Figure  26.  Half-width  at  half-height  for  the  standard  deviation. 
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Figure  27.  An  example  of  a  histogram  for  a  set  of  white  FM  noise. 


5.18.  Frequency  Offset 

It  is  often  necessary  to  estimate  the  frequency  offset  from  either  phase  or  frequency  data. 
Frequency  offset  is  usually  calculated  from  phase  data  by  either  of  three  methods: 
1.         A  least  squares  linear  fit  to  the  phase  data  (optimum  for  white  PM  noise): 
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x(t)  =  a  +  bt,  where  slope  =  y(t)  =  b. 


2.  The  average  of  the  first  differences  of  the  phase  data  (optimum  for  white  FM  noise): 
y(t)  =  slope  =  [x(t+T)  -  x(t)]/T. 

3.  The  difference  between  the  first  and  last  points  of  the  phase  data: 

y(t)  =  slope  =  [  x(end)  -  x(start)  ]  /  (M-1),  where  M  =  #  phase  data  points. 
This  method  is  used  mainly  to  match  the  two  endpoints. 

5.19.  Frequency  Drift 

Most  frequency  sources  have  frequency  drift,  and  it  is  often  necessary  (and  usually  advisable)  to  characterize  and 
remove  this  systematic  drift  before  analyzing  the  stochastic  noise  of  the  source.  The  term  drift  refers  to  the  systematic 
change  in  frequency  due  to  all  effects,  while  aging  includes  only  those  effects  internal  to  the  device.  Frequency  drift 
is  generally  analyzed  by  fitting  the  trend  of  the  frequency  record  to  an  appropriate  mathematical  model  (e.g.,  linear, 
log,  etc.),  often  by  the  method  of  least  squares.  The  model  may  be  based  on  some  physical  basis  or  simply  a 
convenient  equation,  using  either  phase  or  frequency  data,  and  its  suitability  may  be  judged  by  the  degree  to  which  it 
produces  white  (i.e.,  uncorrelated)  residuals. 

Frequency  drift  is  the  systematic  change  in  frequency  due  to  all  effects,  while  frequency  aging  is  the  change  in 
frequency  due  to  effects  within  the  device.  Thus,  for  a  quartz  crystal  oscillator,  aging  refers  to  a  change  in  the 
resonant  frequency  of  its  quartz  crystal  resonator,  while  drift  would  also  include  the  effects  of  its  external 
environment.  Therefore,  drift  is  the  quantity  that  is  usually  measured,  but  it  is  generally  done  under  constant 
environmental  conditions  to  the  greatest  extent  possible  so  as  to  approximate  the  aging  of  the  device. 

5.20.  Drift  Analysis  IVIethods 

Several  drift  methods  are  useful  for  phase  or  frequency  data  as  described  below.  The  best  method  depends  on  the 
quality  of  the  fit,  which  can  be  judged  by  the  randomness  of  the  residuals. 


Table  13.  Drift  analysis  methods  for  phase  or  frequency  data. 


Data 

Method 

Noise  Model 

Phase 

Quadratic  Fit 

W  PM 

Phase 

Avg  of  2nd  Diffs 

RW  FM 

Phase 

3-Point  Fit 

W  &  RW  FM 

Phase 

Linear  Fit 

Frequency  Offset 

Phase 

Avg  of  1st  Diffs 

Frequency  Offset 

Phase 

Endpoints 

Frequency  Offset 

Freq 

Linear  Fit 

W  FM 

Freq 

Bisection  Fit 

W  &  RW  FM 

Freq 

Log  Fit 

Stabilization 

Freq 

Diffusion  Fit 

Diffusion 

5.21.  Phase  Drift  Analysis 

Three  methods  are  commonly  used  to  analyze  frequency  drift  in  phase  data: 
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1. 


A  least  squares  quadratic  fit  to  the  phase  data: 

x(t)  =  a  +  bt  +  cV,  where  y(t)  =  x'(t)  =  b  +  2ct,  slope  =  y'(t)  =  2c. 

This  continuous  model  can  be  expressed  as  x„  ^  a  +  robn  +  vo^cn^  for  n  =1,  2,  3...,  N,  and  to  is  the  sampling 
interval  for  discrete  data,  where  the  a,  b  and  c  coefficients  have  units  of  sec,  sec/sec  and  sec/sec',  respectively, 
and  the  frequency  drift  slope  and  intercept  are  2c  and  b,  respectively  .  The  fit  coefficients  can  be  estimated  by 
the  following  expressions  [1]: 
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where  the  A-F  terms  are  as  follows: 


.4  =  3[3^(A^  +  l)  +  2] 
5  =  -18(2^  +  1) 
C  =  30 

Z)  =  1 2(2^  +  1)(8A^  +  1 1)  /  [(^  +  l)(yV  +  2)] 

£  =  -180/(A^  +  2) 

F  =  180/[(yV  +  l)(A^  +  2)] 

G^N[N  -\)[N  -2) 

A  quadratic  fit  to  the  phase  data  is  the  optimum  model  for  white  PM  noise. 

2.  The  average  of  the  second  differences  of  the  phase  data: 

y(t)  =  [x(t+T)-x(t)]/T,  slope  =  [y(t+T)-y(t)  ]/t  =  [x(t+2T)-2x(t+T)+x(t)]/T^ 
This  method  is  optimum  for  random  walk  FM  noise. 

3.  A  three-point  fit  at  the  start,  middle,  and  end  of  the  phase  data: 

slope  -  4[x(end)-2x(mid)+x(start)]/(MT)',  where  M  =  the  number  of  data  points. 
It  is  the  equivalent  of  the  bisection  method  for  frequency  data. 

5.22.  Frequency  Drift  Analysis 

Four  methods  are  commonly  used  to  analyze  frequency  drift  in  frequency  data: 
1.         A  least  squares  linear  regression  to  the  frequency  data: 
y(t)  =  a  +  bt,  where  a  =  intercept,  b  =  slope  =  y'(t). 
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Linear  frequency  drift  can  be  estimated  by  a  linear  least  squares  fit  to  the  frequency  data,  y(t)  =  a  +  bt.  That 
continuous  model  can  be  expressed  as  v„  =  a  +  mbn  for  n  =/,  2,  3...,  M  where  M  is  the  number  of  frequency 
data  points,  ro  is  the  sampling  interval  for  the  discrete  data.  The  frequency  drift  intercept  and  slope  are  a  and 
b,  and  have  units  of  sec  and  sec/sec,  respectively.  The  fit  coefficients  can  be  estimated  by  equations  55  and 
56: 
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A  linear  fit  to  the  frequency  data  is  the  optimum  model  for  white  FM  noise. 

2.  The  frequency  averages  over  the  first  and  last  halves  of  the  data: 

slope  =  2  [  y(2nd  half)  -  y(lst  half)  ]  /  CNi),  where  N  =  number  of  points. 
This  bisection  method  is  optimum  for  white  and  random  walk  FM  noise. 

3.  A  log  model  of  the  form  (see  MIL-0-553 1  OB)  that  applies  to  frequency  stabilization: 


y(t)  =  a-ln(bt+l ),  where  slope  =  y'(t)  =  ab/(bt+l ). 
4.        A  diffusion  (Vt)  model  of  the  form 

y(t)  =  a+b(t+c)"-,  where  slope  =  y'(t)  =  '/2-b(t+c)""-. 
References  for  Frequency  Drift 


1.  J. A.  Barnes,  "The  Measurement  of  Linear  Frequency  Drift  in  Oscillators,"  Proc.  15th  Anmi.  PTTI  Meeting,  pp. 
551-582,  Dec.  1983. 

2.  J.  A.  Barnes,  "The  Analysis  of  Frequency  and  Time  Data,"  Austron,  Inc.,  Dec.  1991. 

3.  M.A.  Weiss  and  C.  Hackman,  "Confidence  on  the  Three-Point  Estimator  of  Frequency  Drift,"  Proc.  24th  Annu. 
PTTI  Meeting,  pp.  451-460,  Dec.  1992. 


5.23.  All  Tau 


Stability  calculations  made  at  all  possible  tau  values  can  provide  an  excellent  indication  of  the  variations  in  the 
results,  and  are  a  simple  form  of  spectral  analysis.  In  particular,  cyclic  variations  are  often  the  result  of  interference 
between  the  sampling  rate  and  some  periodic  instability  (such  as  environmental  sensitivity).  However,  an  all  tau 
analysis  is  computationally  intensive  and  can  therefore  be  slow.  For  most  purposes,  however,  it  is  not  necessary  to 
calculate  values  at  every  tau,  but  instead  to  do  so  at  enough  points  to  provide  a  nearly  continuous  curve  on  the  display 
device  (screen  or  paper).  Such  a  "many  tau"  analysis  can  be  orders  of  magnitude  faster  and  yet  provide  the  same 
information. 
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(b)  Many  Tau  Stability  Plot 
Figure  28.  Comparison  of  all  tau  and  many  tau  stability. 


5.25.  Environmental  Sensitivity 


Environmental  sensitivity  should  be  treated  separately  from  noise  when  one  performs  a  stability  analysis.  However,  it 
can  be  very  difficult  to  distinguish  between  those  different  mechanisms  for  phase  or  frequency  variations.  It  is  often 
possible  to  control  the  environmental  conditions  sufficiently  well  during  a  stability  run  so  that  environmental  effects 
such  as  temperature  and  supply  voltage  are  negligible.  Determining  how  well  those  factors  have  to  be  controlled 
requires  knowledge  of  the  device's  environmental  sensitivities.  Those  factors  should  be  measured  individually,  if 
possible,  over  the  largest  reasonable  excursions  to  minimize  the  effect  of  noise.  Environmental  sensitivity  can  best  be 
determined  by  considering  the  physical  mechanisms  that  apply  within  the  unit  under  test.  Useful  information  about 
the  environmental  sensitivity  of  frequency  sources  can  be  found  in  the  references  below.  Some  environmental  factors 
affect  phase  and  frequency  differently,  which  can  cause  confusion.  For  example,  temperature  affects  the  phase  delay 
through  a  cable.  Dynamically,  however,  a  temperature  ramp  produces  a  rate  of  change  of  phase  that  mimics  a 
frequency  change  in  the  source.  Because  environmental  sensitivity  is  highly  dependent  on  device  and  application,  it 
does  not  receive  detailed  consideration  in  this  handbook.  More  information  will  be  found  in  the  following  references. 
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References  for  Environmental  Sensitivity 

1.  IEEE  Standard  1 193.  IEEE  Guide  for  Measurement  of  Environmental  Sensitivities  of  Standard  Frequency 
Generators. 

2.  MIL-0-55310.  General  Specification  for  Military  Specification,  Oscillators,  Crystal,  Military  Specifications 
and  Standards. 

3.  MIL-STD-8 1 0,  Environmental  Test  Methods  and  Engineering  Guidelines. 

4.  MIL-STD-202,  Test  Methods  for  Electronic  and  Electrical  Component  Parts. 

5.  C.  Audion,  N.  Dimarcq,  V.  Giodano,  J.  Viennet,  "Physical  Origin  of  the  Frequency  Shifts  in  Cesium  Beam 
Frequency  Standards:  Related  Environmental  Sensitivity,"  Proceedings  of  the  22'"'  Annual  Precise  Time  and 
Time  Inter\-al  (PTTI)  Applications  and  Planning  Meeting,  pp.  419-440,  1990. 

6.  L.A.  Breakiron.  L.S.  Cutler.  H.W.  Hellvvig,  J.R.  Vig,  G.M.R  Winkler,  "General  Considerations  in  the 
Metrology  of  the  Environmental  Sensitivities  of  Standard  Frequency  Generators,"  Proceedings  of  the  1992 
IEEE  Frequency  Control  Symposium,  pp.  816-830,  1992. 

7.  H.  Hellwig,  "Environmental  Sensitivities  of  Precision  Frequency  Sources,"  IEEE  Transactions  lM-39,  pp. 
301-306.  r990. 

8.  E.M.  Mattison,  "Physics  of  Systematic  Frequency  Variations  in  Hydrogen  Masers,"  Proceedings  of  the  22"'' 
Annual  Precise  Time  and  Time  Interval  (PTTI)  Applications  and  Planning  Meeting,  pp.  453-464.  1990. 

9.  W.J.  Riley.  "The  Physics  of  Environmental  Sensitivity  of  Rubidium  Gas  Cell  Atomic  Frequency  Standards." 
Proceedings  of  the  22"'^  Annual  Precise  Time  and  Time  Inter\'al  (PTTI)  Applications  and  Planning  Meeting, 
pp.  441-452,  1990. 

10.  F.L.  Walls  and  J.J.  Gagnepain,  "Environmental  Sensitivities  of  Quartz  Crystal  Oscillators,"  IEEE 
Transactions  UFFC-39,  no.  2,  pp.  241-249,  1992. 

11.  P.S.  Jorgensen.  "Special  Relativity  and  Intersatellite  Tracking,"  Ncnngation,  Vol.  35,  No.  4.  pp.  429-442. 
Winter  1^988-89. 

12.  T.E.  Parker,  "Environmental  Factors  and  Hydrogen  Maser  Frequency  Sensitivity,"  IEEE  Trans.  UFFC-46, 
No.  3.  pp.745-751.  1999. 


5.26.  Parsimony 

In  any  measurement  or  analysis,  it  is  desirable  to  minimize  the  number  of  extraneous  parameters.  This  is  not  just  a 
matter  of  elegance;  the  additional  parameters  may  be  implicit  or  arbitrary  and  thereby  cause  confusion  or  error  if  they 
are  ignored  or  misunderstood.  For  example,  the  Allan  deviation  has  the  important  advantage  that,  because  of  its 
convergence  characteristics  for  all  common  clock  noises,  its  expected  value  is  independent  of  the  number  of  data 
points.  Many  of  the  techniques  used  in  the  analysis  of  frequency  stability,  however,  do  require  that  certain  parameters 
be  chosen  before  they  can  be  performed.  For  example,  drift  removal  requires  the  choice  of  a  mode!  to  which  the  data 
will  be  fit  (perhaps  using  the  criterion  of  white  residuals).  Outlier  removal  is  an  especially  difficult  case,  where 
judgment  often  enters  into  the  decision  as  to  whether  or  not  a  datum  is  anomalous.  A  listing  of  some  of  these 
nonparsimonious  parameters  is  given  in  Table  14: 
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Table  14.  Non-Parsimonious  Parameters  in  Frequency  Stability  Analysis 


Type 

Process 

Parameter 

Criterion 

Remarks 

Pre- 

Outlier 

Number  of 

Apply  best 

Use  of  MAD-based 

processing 

removal 

sigmas 

judgment 

robust  statistics  is 
recommended. 

Drift 

Remove  or 

Is  drift 

It  is  generally  wise 

removal 

not 

deterministic? 
Is  its  cause 
known? 

to  remove 
deterministic  drift 
beiore  noise 
analysis 

Model 

White 

residuals  are 
sign  that  model 
is  appropriate. 

Model  may  have 
physical  basis  (e.g., 
ditrusion  process) 

Convergence 

For  iterative 

Generally  uncritical 

limit 

tit. 

-  can  hide  deeply  in 
algorithm 

Remove 

Model 

Noise  type 

Not  necessarily 

average 

simple  arithmetic 

frequency 

mean  of  frequency. 

Phase  to 

Tau 

Accuracy 

Is  measurement 

frequency 

interval  known 

conversion 

exactly?  Is  it  the 
same  for  each  point? 

Frequency 

Tau 

See  above 

to  phase 

Initial  phase 

Set  to  zero. 

Generally  arbitrary. 

conversion 

Doesn't  affect 
subsequent  stability 
analysis. 

Normalize 

Use  to 

See  "'Remove 

frequency 

emphasize 
noise  rather 
than  frequency 
offset 

average  frequency 
above" 

Analysis 

Drift 

Model 

Smallest 

Can  be  critical. 

estimation 

residuals 

especially  for 
predictions.  Known 
physics  can  help 
choose. 

Convergence 

For  iterative  fit 

Generally  uncritical 

limit 

-  can  hide  deeply  in 
algorithm 

Frequency 

Noise  model 

Lowest 

Generally  uncritical 

estimation 

residuals  for 
known  noise 
type 

Gaps 

Skip,  fill, 

Number, 

Process  choice 

omit,  or 

distribu-tion. 

affects  results 

exclude 

noise  type 
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Allan 

Number  of 

Time  available 

As  many  as  possible 

deviation 

data  points 

tor 

measurement. 
Must  remove 
outliers. 

Dead  time 

Property  of 
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5.27.  Transfer  Functions 

Variances  can  be  related  to  the  spectral  density  of  the  fractional  frequency  fluctuations  by  their  transfer  functions.  For 
example,  for  the  Hadamard  variance,  this  relationship  is 

^\{r)=l'S^{f)\H,{f)Uf^  (57) 

where  <5^w{x)  is  the  three-sample,  zero  dead-time,  binominally  weighted  Hadamard  variance,  Sy(/)  is  the  spectral 
density  of  the  fractional  frequency  fluctuations,  Hh(/)  is  its  transfer  function,  and  /h  is  the  upper  cutoff  frequency 
(determined  by  hardware  factors).  The  transfer  function  is  determined  by  the  statistic's  time-domain  sampling 
pattern. 

The  transfer  functions  for  the  most  common  variances  used  in  frequency  stability  analysis  are  shown  in  Table  15: 


Table  15.  Transfer  functions  for  the  most  common  variances. 


Variance 

Magnitude  Squared  Transfer  Function  |H(f)|^ 

Allan 

2 

'  ■  ^\ 
sin  KTJ 

.  ) 

2 

sin"  Kzf 

Hadamard 

.  2' 

(  ■ 
sin  KTJ 

sin'*  TTzf 
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For  nxf «  1,  the  transfer  function  of  the  Allan  variance  behaves  as  {nxfy,  indicating  that  it  is  convergent  for  power 
law  processes  Sy"  down  to  as  low  as  a  =  -2  (Random  Walk  FM),  while  the  transfer  function  of  the  Hadamard  variance 
behaves  as  (ttt/)^,  indicating  that  it  is  convergent  for  power  law  processes  as  low  as  a  =  -4  (Random  Run  FM). 

The  squared  magnitudes  of  these  transfer  functions  are  shown  in  the  plots  below: 


Transfer  Function  of  Allan  Variance 
2-Sample  (N=2),  Zero  Dead-Time  (T=i) 


— <  0.6 


(a)  Allan  Variance 


Transfer  Function  of  Hadamard  Vanance 
3-Sample  {N=3).  Zero  Dead-Time  (T=t),  Binomially  Weighted  Coefficients  (1,2,1) 


|H„(f) 

'=2''(sin{ 

TTf)/Rtf)'( 

5in''(itTf) 

(b)  Hadamard  Variance 

Figure  29.  Squared  magnitudes  of  transfer  functions  for  (a)  Allan  Variance  and  (b)  Hadamard  Variance. 

These  responses  have  their  peaks  where  the  frequency  is  one-half  the  sampling  rate,  and  nulls  where  it  is  a  multiple  of 
the  sampling  rate  (i.e.,  at  f  =  n/i,  where  n  is  an  integer).  As  a  spectral  estimator,  the  Hadamard  variance  has  slightly 
higher  resolution  than  the  Allan  variance,  since  the  equivalent  noise  bandwidths  of  the  Hadamard  and  Allan  spectral 
windows  are  0.41 1  i''  and  0.476  t  '.  respectively  [5]. 


Similar  transfer  functions  exist  for  the  modified,  total,  and  Theol  variances. 
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6 


Frequency  Domain  Stability 


Frequency  stability  can  also  be  characterized  in  the  frequency  domain  in  terms 
of  a  power  spectral  density  (PSD)  that  describes  the  intensity  of  the  phase  or 
frequency  fluctuations  as  a  function  of  Fourier  frequency.  Spectral  stability 
measures  are  directly  related  to  the  underlying  noise  processes,  and  are 
particularly  appropriate  when  the  phase  noise  of  the  source  is  of  interest. 

6.1.   Noise  Spectra 


Frequency  domain  stability 
measures  are  based  on  power 
spectral  densities  that  characterize 
the  intensity  of  the  phase  or 
frequency  fluctuations  as  a 
function  of  Fourier  frequency. 


The  random  phase  and  frequency  fluctuations  of  a  frequency  source  can  be  modeled  by  power  law  spectral  densities 
of  the  form 


Sv(f)  =  h(a)f\ 

where;  Sy(f)     =         one-sided  power  spectral  density  of  y,  the 
fractional  frequency  fluctuations,  1/Hz 
f  =         Fourier  or  sideband  frequency,  hertz, 

h(a)      =  intensity  coefficient,  and 

a         =         exponent  of  the  power  law  noise  process. 


The  most  commonly  encountered  noise  spectra  are 


White  (f") 
Flicker  (f-') 
Random  Walk(f--) 
Flicker  Walk  (f  ^) . 


Examples  of  these  noise  types  are  shown  in  the  figure  below. 
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Figure  30.  Examples  of  noise  types. 

Power  law  spectral  models  can  be  applied  to  both  phase  and  frequency  power  spectral  densities.  Phase  is  the  time 
integral  of  frequency,  so  the  relationship  between  them  varies  as  1/f : 
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S.(f)-7T^,  (59) 
where   Sx(f)     =         PSD  of  the  time  fluctuations,  sec^/Hz.  (60) 

Two  other  quantities  are  also  commonly  used  to  measure  phase  noise: 

S^{f)  =         PSD  of  the  phase  fluctuations,  radVHz  and  its 

logarithmic  equivalent  £(f),  dBc/Hz. 


The  relationship  between  these  is 


and 


(61) 


a/)  =  10 -log 


■s,(/) 


where  vq  is  the  carrier  frequency,  hertz. 


(62) 


The  power  law  exponent  of  the  phase  noise  power  spectral  densities  is  p=  a-2.  These  frequency-domain  power  law 
exponents  are  also  related  to  the  slopes  of  the  following  time  domain  stability  measures: 


Allan  variance  cr-y(T)  n=  -(a+1),  a<2 

Modified  Allan  variance  Mod  cry(i)  =  -(a+l,a<3 

Time  variance  <^~\{'^)  n  ^  -(a-1),  a<3. 


The  spectral  characteristics  of  the  power  law  noise  processes  commonly  used  to  describe  the  performance  of 
frequency  sources  are  shown  in  the  following  table: 


Table  16.  Spectral  characteristics  of  po 


wer  law  noise  processes 


Noise  Type 

a 

P 

^' 

White  PM 

2 

0 

-2 

-3 

-1 

Flicker  PM 

1 

-1 

-2 

-2 

0 

White  FM 

0 

_2 

-1 

-1 

1 

Flicker  FM 

-1 

-3 

0 

0 

2 

Random  walk  FM 

-2 

-4 

1 

1 

3 

6.2.   Power  Spectral  Densities 


Four  types  of  power  spectral  density  are  commonly  used  to  describe  the  stability  of  a  frequency  source: 
PSD  of  Frequency  Fluctuations  Sy(0 

The  power  spectral  density  of  the  fractional  frequency  fluctuations  y(t)  in  units  of  1/Hz  is  given  by  Sy(f)  =  h(a)  •  f  °, 
where  f=  sideband  frequency,  Hz. 


PSD  of  Phase  Fluctuations  S^O 
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The  power  spectral  density  of  the  phase  fluctuations  in  units  of  radVHz  is  given  by  S^(f) 
carrier  frequency,  Hz. 


=  (27rv,i)'  •  Sx(f),  where  v„  = 


PSD  of  Time  Fluctuations  S,(f) 

The  power  spectral  density  of  the  time  fluctuations  x(t)  in  units  of  secVHz  is  given  by  Sx(0  =  h(p)  ■  f "  =  Sy{0/(27tf)', 
where:  p=  a-2.  The  time  fluctuations  are  related  to  the  phase  fluctuations  by  x(t)  =  (|)(t)/2rtVo.  Both  are  commonly 
called  "phase"  to  distinguish  them  from  the  independent  time  variable,  t. 

SSB  Phase  Noise  £(f) 

The  SSB  phase  noise  in  units  of  dBc/Hz  is  given  by  £(f)  =  10  •  log  ['A  ■  S^(f)].  This  is  the  most  common  function  used 
to  specify  phase  noise. 

6.3.   Phase  Noise  Plot 

The  following  diagram  shows  the  slope  of  the  SSB  phase  noise,  £(f),  dBc/Hz  versus  log  f  Fourier  frequency,  Hz  for 
various  power  law  noise  processes. 


SSB  Phase  Noise  Diagram 
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Figure  31.  SSB  Phase  Noise  Plot. 


6.4.   Spectral  Analysis 

Spectral  analysis  is  the  process  of  characterizing  the  properties  of  a  signal  in  the  frequency  domain,  either  as  a  power 
spectral  density  for  noise,  or  as  the  amplitude  and  phase  at  discrete  frequencies.  Spectral  analysis  can  thus  be  applied 
to  both  noise  and  discrete  components  for  frequency  stability  analysis.  For  the  former,  spectral  analysis  complements 
statistical  analysis  in  the  time  domain.  For  the  latter,  spectral  analysis  can  aid  in  the  identification  of  periodic 
components  such  as  interference  and  environmental  sensitivity.  Time  domain  data  can  be  used  to  perform  spectral 
analysis  via  the  Fast  Fourier  Transform  (FFT),  and  there  is  much  technical  literature  on  that  subject  [2,3].  While,  in 
principle,  time  and  frequency  domain  analyses  provide  equivalent  information,  in  practice,  one  or  the  other  is  usually 
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preferred,  either  because  of  measurement  and/or  analysis  convenience,  or  because  the  results  are  more  applicable  to  a 
particular  application.  Spectral  analysis  is  most  often  used  to  characterize  the  short-term  (<  1  s)  fluctuations  of  a 
frequency  source  as  a  plot  of  phase  noise  power  spectral  density  versus  sideband  frequency,  while  a  time  domain 
analysis  is  most  often  used  to  provide  information  about  the  statistics  of  its  instability  over  longer  intervals  (>  1  s). 
Modern  instrumentation  is  tending  to  merge  these  approaches  by  digitizing  the  signal  waveform  at  high  sampling 
rates,  thereby  allowing  FFT  analysis  at  relatively  high  Fourier  frequencies.  Nevertheless,  there  are  many  pitfalls  and 
subtleties  associated  with  spectral  analysis  that  must  be  considered  if  meaningful  results  are  to  be  obtained. 


Data  windowing  is  the  process  of  applying  a  weighting  function  that  falls  off  smoothly  at  the  beginning  and  end  to 
avoid  spectral  leakage  in  an  FFT  analysis.  Without  windowing,  bias  will  be  introduced  that  can  severely  restrict  the 
dynamic  range  of  the  PSD  result.  The  most  common  windowing  types  are  Manning,  Hamming,  and  Multitaper.  The 
classic  Manning  and  Hamming  windows  can  be  applied  more  than  one  time. 


Without  filtering  or  averaging,  the  variances  of  the  PSD  results  are  always  equal  to  their  values  regardless  of  the  size 
of  the  time  domain  data  set.  More  data  provide  finer  frequency  resolution,  not  lower  noise  (while  the  data  sampling 
time  determines  the  highest  Fourier  frequency).  Without  averaging,  for  white  noise,  each  spectral  result  has  only  two 
degrees  of  freedom.  Some  sort  of  filtering  or  averaging  is  usually  necessary  to  provide  less  noise  in  the  PSD  results. 
This  can  be  accomplished  by  dividing  the  data  into  sections,  performing  an  FFT  analysis  on  each  section  separately, 
and  then  averaging  those  values  to  obtain  the  final  PSD  result.  The  averaging  factor  improves  the  PSD  standard 
deviation  by  the  square  root  of  the  averaging  factor.  The  tradeoff  in  this  averaging  process  is  that  each  section  of  the 
data  is  shorter,  yielding  a  result  with  coarser  frequency  resolution  that  does  not  extend  to  as  low  a  Fourier  frequency. 


The  multitaper  PSD  analysis  method  offers  a  better  compromise  among  bias,  variance  and  spectral  resolution. 
Averaging  is  accomplished  by  applying  a  set  of  orthogonal  windowing  (tapering)  functions  called  discrete  prolate 
spheroidal  sequences  (DPSS)  or  Slepian  functions  to  the  entire  data  array.  An  example  of  seven  of  these  functions  for 
order  J=4  is  shown  in  Figure  32. 


6.5.   PSD  Windowing 


6.6.   PSD  Averaging 


6.7.   IVIultitaper  PSD  Analysis 
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Figure  32.  Slepian  SPSS  taper  functions,  J=4,  #=7. 
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The  first  function  resembles  a  classic  window  function,  while  the  others  sample  other  portions  of  the  data.  The  higher 
windows  have  larger  amplitude  at  the  ends  that  compensates  for  the  denser  sampling  at  the  center.  These  multiple 
tapering  functions  are  defined  by  two  parameters,  the  order  of  the  function,  J,  which  affects  the  resolution  bandwidth, 
and  the  number  of  windows,  which  affects  the  variance.  A  higher  J  permits  the  use  of  more  windows  without 
introducing  bias,  which  provides  more  averaging  (lower  variance)  at  the  expense  of  lower  spectral  resolution,  as 
shown  Table  17; 


Table  17.  The  two  parameters  of  the  tapering  function. 


Order  J 

#  Windows 

2.0 

1-3 

2.5 

1-4 

3.0 

1-5 

3.5 

1-6 

4.0 

1-7 

4.5 

1-8 

5.0 

1-9 

The  resolution  BW  is  given  by  2J/Nt,  where  N  is  the  number  of  data  points  sampled  at  time  interval  /.  An  adaptive 
algorithm  can  be  used  to  weight  the  contributions  of  the  individual  tapers  for  lowest  bias.  The  multitaper  PSD  has  a 
flat-topped  response  for  discrete  spectral  components  that  is  nevertheless  narrower  than  an  averaged  periodogram 
with  the  same  variance.  It  is  therefore  particularly  useful  for  examining  discrete  components  along  with  noise. 

6.8.    PSD  Notes 

A  carrier  frequency  parameter  applies  to  the  S^iV)  and  £(f)  PSD  types.  The  number  of  Fourier  frequency  points  is 
always  the  power  of  2  greater  than  or  equal  to  one-half  of  the  number  of  time  domain  data  points,  n.  The  spacing 
between  Fourier  frequency  points  is  1/nt,  and  the  highest  Fourier  frequency  is  l/2t.  If  averaging  is  done,  the  value  of 
n  is  reduced  by  the  averaging  factor.  The  PSD  fit  is  a  least-squares  power  law  line  through  octave-band  PSD  averages 
[6]. 

For  characterizing  frequency  stability,  a  spectral  analysis  is  used  primarily  for  the  analysis  of  noise  (not  discrete 
components),  and  should  include  the  quantitative  display  of  power  law  noise  in  common  PSD  units,  perhaps  with  fits 
to  integer  power  law  noise  processes.  Amplitude  corrections  need  to  be  made  for  the  noise  response  of  the 
windowing  functions.  The  amplitude  of  discrete  components  should  be  increased  by  the  log  of  the  BW  (Fourier 
frequency  spacing  in  hertz),  which  is  a  negative  number  for  typical  sub-hertz  bandwidths. 
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7      Domain  Conversions 


The  stabilit\  of  a  frequency  source  can  be  specified  and  measured  in  either  the 
time  or  the  frequency  domain.  Examples  of  these  stability  measures  are  the 
Allan  variance.  a\(T),  in  the  time  domain,  and  the  spectral  densitv  of  the 
fractional  frequency  fluctuations,  Sy(f),  in  the  frequency  domain.  Conversions 
between  these  domains  may  be  made  by  numerical  integration  of  their 
fundamental  relationship,  or  by  an  approximation  method  based  on  a  power  law 
spectral  model  for  the  noise  processes  involved.  The  latter  method  can  be 
applied  onl\  when  the  dominant  noise  process  decreases  toward  higher 
sideband  frequencies.  Otherwise,  the  more  fundamental  method  based  on 
numerical  integration  must  be  used.  The  general  conversion  from  time  to 
frequency  domain  is  not  unique  because  white  and  flicker  phase  noise  have  the 
same  Allan  variance  dependence  on  t.  When  performing  any  of  these 
conversions,  it  is  necessar\  to  choose  a  reasonable  range  for  a  and  i  in  the 
domain  being  converted  to.  The  main  lobe  of  the  ay(i)  and  Mod  UyCx)  responses 
occur  at  the  Fourier  frequenc>  f  =  1/2t. 

Time  domain  frequency  stabilit}  is  related  to  the  spectral  density  of  the  fractional  frequencv  fluctuations  b\  the 
relationship 


The  stability  of  a  frequency  source 
can  be  specified  and  measured  in 
either  the  time  or  the  frequency 
domain.  One  domain  is  often 
preferred  to  specify  the  stability 
because  it  is  most  closely  related 
to  the  particular  application. 
Similarly,  stability  measurements 
may  be  easier  to  perform  in  one 
domain  than  the  other. 
Conversions  are  possible  between 
these  generally  equivalent 
measures  of  frequency  stability. 


ctHM.T.t)^  ^'S^if)-\Hif)f-df, 


(63) 


where  ^H(fj  '  is  the  transfer  function  of  the  time  domain  sampling  function. 

The  transfer  function  of  the  Allan  (two-sample)  time  domain  stability  is  given  by 


\H{ft=2 


(64) 


and  therefore  the  Allan  variance  can  be  found  from  the  frequency  domain  by  the  expression 


a]{T)  =  2[s^{f) 


The  equivalent  expression  for  the  modified  Allan  variance  is 


(65) 


Mod  (7 


2  ^f.SSf)s\n\7rTf)^^ 

;(^)=  ,-4  -  ^  I      .  V — TT^J 


(66) 


7.1.   Power  Law  Domain  Conversions 

Domain  conversions  may  be  made  for  power  law  noise  models  by  using  the  following  equation  and  conversion 
formulae: 


2.  ^    u  ^01      o    ^,   ^   .     1 .038  +  3  log^.  2;r./;,r  ,  ,  3/, 

^6  2z  2k  T  \2k)  r 


(67) 


73 


where  the  h„  terms  define  the  level  of  the  various  power  law  noises. 


Noise  Type  g^vd) 


RW  FM  a'  f  -  •  Sv(f)  ■  T  ' 

F  FM  B  •  f  '  •  S,(0  ■i'^ 

W  FM  C  •  f°  •  Sy{f)  •  T-' 

FPM  D  •  f  -' •  Sv(0  ■  T  ^ 

WPM  E-f-Sy(tVr- 


-t"-  aMt)-f-' 
C-'  -i'  •  ayT)-f° 
D-'  •  I  -  ■  aMi)  -  f  ' 
E    ■  X  '  ■  a^t)  •  f 


where 


A  =  47tV6 
B  =  21n2 


C  =  1/2 

D=  1.038  +  3-ln(2irfhTo)/47i^ 
E  =  3fh/47:^ 

and  fh  is  the  upper  cutoff  frequency  of  the  measuring  system  in  hertz,  and  tq  is  the  basic  measurement  time.  This 
factor  applies  only  to  white  and  flicker  PM  noise.  The  above  conversion  formulae  apply  to  the  TheoH  hybrid  statistic 
as  well  as  to  the  Allan  variance. 

7.2.   Example  of  Domain  Conversions 

This  section  shows  an  e.xample  of  time  and  frequency  domain  conversions.  First,  a  set  of  simulated  power  law  noise 
data  is  generated,  and  the  time  domain  properties  of  this  noise  are  analyzed  by  use  of  the  overlapping  Allan  deviation. 
Next,  the  same  data  are  analyzed  in  the  frequency  domain  with  an  £(f)  PSD.  Then,  a  power  law  domain  conversion  is 
done,  and  those  results  are  compared  with  those  of  the  spectral  analysis.  Finally,  the  other  power  spectral  density 
types  are  examined. 

For  this  example,  we  generate  4097  phase  data  points  of  simulated  white  FM  noise  with  a  one-second  Allan  deviation 
value  ay(l)  =  1  x  10"  and  a  sampling  interval  t  =  1  ms.  The  number  of  points  is  chosen  as  an  even  power  of  2  for 
the  frequency  data.  The  generated  set  of  simulated  white  FM  noise  is  shown  as  frequency  data  in  Figure  33a,  and  their 
overlapping  Allan  deviation  is  shown  in  Figure  33b.  The  ay(l)white  FM  noise  fit  parameter  is  1.08e-l  1,  close  to  the 
desired  value. 
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FREQUENCY  DATA 

Simj^ted  Whits  FU  None 


2.0  3.0 

Time,  Seconds 


FREQUENCY  STABILITY 

SlmulatBd  Wtilts  FM  Haiso 


1.00^-03  3.1S!:-10 
2.00e-0J  ;.23e-lO 
•l.OOc-OJ  1.58C-I0 
8.00f-05  1.16<-ig 
1.60e-02  8.0I--11 
3,20p-02  6.(X)e-lI 
6.'0e-02  4.99e-l1 
l.28e-CH  3.35e-11 
2.56«-01  1.95T-11 
5.t2e-01  1.47C-11 
1.025+00  1.23^-11 
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Av«a|^^  Time,  Seconds 


10" 


(a) 


(b) 


POWER  SPECTRUM 

Simulate  Whits  FM  No>>« 


"To"  io"  W~ 

Fourier  Fieqoeiicy,  f.  Hz 


-6 

N  -7 

-  -9 

c5 

"  -u 

ao 

2-n 


o 

Z  -14 

a.  -16 


POWER  SPECTRUM 

Simulated  WWta  FM  NoM 


IWndowM  Honnina,  Cofiiar  frjq=:1.00O«-t-007[ 


lO"  Iff  lo' 

Fourier  Frequency,  f,  Hz 


105 


(C] 


POWER  SPECTRUM 

SImuJated  White  FM  NoKe 


POWER  SPECTRUM 

simulated  Whits  FM  Noise 


Fourier  Frequency,  f,  Hz 


Fourier  Frequency,  f,  Hz 


(e) 


(f) 


Figure  33.  (a)  Simulated  W  FM  Noise  Data,  (b)  Overlapping  Allan  Deviation,  (c)  £(f)  Power  Spectral  Density,  (d)  S^(f) 
Power  Spectral  Density,  (e)  Sx(f)  Power  Spectral  Density,  (f)  Sy(f)  Power  Spectral  Density 


The  power  spectrum  for  the  phase  data  is  calculated  by  use  of  a  10  MHz  carrier  frequency  and  a  £(f)  power  spectral 
density  type,  the  SSB  phase  noise  to  signal  ratio  in  a  1  Hz  BW  as  a  function  of  sideband  frequency,  f,  as  shown  in 
Figure  33c.  The  fit  parameters  show  an  £(1)  value  of -79.2  dBc/Hz  and  a  slope  of -19.6  dB/decade.  in  close 
agreement  w  ith  the  expected  values  of  -80  dBc/Hz  and  -20  dB/decade. 
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The  expected  PSD  values  that  correspond  to  the  time  domain  noise  parameters  used  to  generate  the  simulated  power- 
law  noise  can  be  determined  by  the  power  law  domain  conversion  formulae  of  Section  7.1,  as  shown  in  Table  18. 

Table  18.  Domain  calculation  results 

Frequency  Domain:  PSD  Type:   L(f),  dBc/Hz 

SB  Freq  (Hz):  l.OOOOOe+00 
Carrier (MHz) :  l.OOOOOe+01 

Time  Domain:  Sigma  Type:  Normal 

Tau   (Sec)  :    1 .  OOOOOe-03 
Avg  Factor:  1000 


Power 

Law  Noise: 

Type 

dB/dec 

PSD 

Type 

Mu 

Sigma 

RWFM 

-40 

None 

RWFM 

+  1 

O.OOOOOe+00 

FFM 

-30 

None 

FFM 

0 

O.OOOOOe+00 

WFM 

-20 

-80.0 

WFM 

-1 

l.OOOOOe-11 

FPM 

-10 

None 

FPM 

-2 

O.OOOOOe+00 

WPM 

0 

None 

WPM 

-2 

O.OOOOOe+00 

All 

-80  .  0 

All 

1 . OOOOOe-11 

The  other  types  of  PSD  that  are  commonly  used  for  the  analysis  of  frequency  domain  stability  analysis  are  S^(f),  the 
spectral  density  of  the  phase  fluctuations  in  radVHz;  Sx(f),  the  spectral  density  of  the  time  fluctuations  in  sec'/Hz;  and 
Sv(f)-.  the  spectral  density  of  the  fractional  frequency  fluctuations  in  units  of  1/Hz.  The  expected  value  of  all  these 
quantities  for  the  simulated  white  FM  noise  parameters  with  ay(l)  =  l.OOe-11,  t  =  l.OOe-3,  and  fo  =  10  MHz  are 
shown  in  the  following  table. 


Table  19.  Expected  value  of  PSD  types  for  the  simulated  white  FM  noise  parameters. 


Parameter 

PSD  Type 

£(f),  dBc/Hz 

S4,(0,  radVHz 

S,(f),  sec'/Hz 

Sy({),  1/HZ 

Data  Type 

Phase 

Phase 

Phase 

Frequency 

Simulated  Value 

-80 

2e-8 

5.066e-24 

2e-22 

Log  Value 

same 

-7.70 

-23.3 

-21.7 

Slope,  dB/decade 

-20 

-20 

-20 

0 

Fit  Value 

-79.2 

1.19e-8 

3.03e-24 

1.40e-22 

Fit  Exponent 

-1.96 

-1.96 

-1.96 

-0.04 
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8      Noise  Simulation 


It  is  valuable  to  have  a  means  of  generating  simulated  power  law  clock  noise  having  the  desired  noise  type  (white 
phase,  flicker  phase,  white  frequency,  flicker  frequency,  and  random  walk  frequency  noise),  Allan  deviation, 
frequency  offset,  frequency  drift,  and  perhaps  a  sinusoidal  component.  This  can  serve  as  both  a  simulation  tool  and  as 
a  way  to  validate  stability  analysis  software,  particularly  for  checking  numerical  precision,  noise  recognition,  and 
modeling.  A  good  method  for  power-law  noise  generation  is  described  in  Reference  8.  The  noise  type  and  time 
series  of  a  set  of  simulated  phase  data  are  shown  in  Table  20: 

Table  20.  Noise  type  and  time  series  for  a  set  of  simulated  phase  data. 


Noise  Type 

Phase  Data  Plot 

Random  walk  FM 

a  =  -2 
Random  run  noise 

j 

/ 

Flicker  FM 
a  =  -l 
Flicker  walk  noise 

White  FM 
a=0 
Random  walk  noise 

Flicker  PM 

a=  1 
Flicker  noise 

iWiii 

1 

r 

1  1' 

1 

White  PM 

a  =  2 

1  1 1 

 1 

1 

White  noise 
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8.1.  White  Noise  Generation 


White  noise  generation  is  straightforward.  One  popular  technique  is  to  first  generate  two  independent  uniformly 
distributed  random  sequences  [1],  and  combine  them  using  the  Box-Muller  transform  [2,3]  to  produce  a  white 
spectrum  with  Gaussian  deviates.  Another  method  is  to  generate  12  independent  random  sequences  uniformly 
distributed  between  0  and  1,  add  them,  and  subtract  6  [4].  This  will,  via  the  central  limit  theorem,  produce  a  Gaussian 
distribution  having  zero  mean  and  unit  variance.  White  noise  can  be  numerically  integrated  and  differenced  to 
transform  it  by  1/f"  and  f^,  respectively,  to  produce  simulated  noise  having  any  even  power  law  exponent. 

8.2.  Flicker  Noise  Generation 

Flicker  noise  is  more  difficult  to  generate  because  it  cannot  be  described  exactly  by  a  rational  transfer  function,  and 
much  effort  has  been  devoted  to  generating  it  [5-9].  The  most  common  methods  involve  linear  filtering  by  RC  ladder 
networks  [5],  or  by  FFT  transformation  [7,9].  The  FFT  method  can  produce  noise  having  any  integer  power  law 
exponent  from  a  =  -2  (RW  FM)  to  a  =  +2  (W  PM)  [7,  8]. 

8.3.  Flicker  Walk  and  Random  Run  Noise  Generation 

The  more  divergent  flicker  walk  FM  (a  =  -3)  and  random  run  FM  (a  =  -4)  power  law  noise  types  may  be  generated 
by  using  the  1/f^  spectral  property  of  a  frequency  to  phase  conversion.  For  example,  to  generate  RR  FM  noise,  first 
generate  a  set  of  RW  FM  phase  data.  Then  treat  this  RW  FM  phase  data  as  frequency  data,  and  convert  it  to  a  new  set 
of  RR  FM  phase  data. 

8.4.  Frequency  Offset,  Drift,  and  Sinusoidal  Components 

Beside  the  generation  of  the  desired  power  law  noise,  it  is  desirable  to  include  selectable  amounts  of  frequency  offset, 
frequency  drift,  and  a  sinusoidal  component  in  the  simulated  clock  data. 
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9      Measuring  Systems 


A  frequency  measuring  system 
with  adequate  resolution  and  a 
low  noise  floor  is  necessary  to 
make  precision  clock 
measurements. 


Frequency  measuring  systems  are  instruments  that  accept  two  or  more  inputs,  one 
of  which  may  be  considered  to  be  the  reference,  and  compare  their  relative  phase  or 
frequencies.  These  systems  can  take  many  forms,  from  the  direct  use  of  a 
frequency  counter  to  elaborate  low-noise,  high-resolution  multichannel  clock 
measuring  systems  with  associated  archival  databases.  They  can  be  custom  built  or 
bought  from  several  organizations  specializing  in  such  systems.     The  most 

important  attribute  of  a  frequency  measuring  system  is  its  resolution,  which,  for  high  performance  devices,  requires  1 
ps/s  (pplO  ")  or  better  resolution,  and  more  elaborate  hardware  than  a  counter.  The  resolution  of  a  digital  frequency, 
period,  or  time  interval  counter  is  determined  mainly  by  its  speed  (clock  rate)  and  the  performance  of  its  analog 
interpolator  (if  any).  That  resolution  generally  improves  linearly  with  the  averaging  time  of  the  measurement,  and  it 
can  be  enhanced  by  preceding  it  w  ith  a  mixer  that  improves  the  resolution  by  the  heterodyne  factor,  the  ratio  of  the  RF 
input  to  the  IF  beat  frequencies.  Noise  is  another  important  consideration  for  a  high-performance  measuring  system 
whose  useful  resolution  may  be  limited  by  its  noise  floor,  the  scatter  in  the  data  when  the  two  inputs  are  driven 
coherently  by  the  same  source.  The  performance  of  the  measuring  system  also  depends  on  the  stability  of  its 
reference  source.  A  low  noise  ovenized  quartz  crystal  oscillator  may  be  the  best  choice  for  a  reference  in  the  short 
term  (1  to  100  s),  while  a  active  hydrogen  maser  generally  provides  excellent  stability  at  averaging  times  out  to 
several  days,  and  cesium  beam  tube  devices  at  longer  averaging  times. 


Three  methods  are  commonly  used  for  making  precise  time  and  frequency  measurements,  as  described  below. 

9.1.   Time  Interval  Counter  Method 

The  time  interval  counter  method  divides  the  two  sources  being  compared  down  to  a  much  lower  frequency  (typically 
1  pulse/second)  and  measures  their  time  difference  with  a  high  resolution  time  interval  counter: 


Comp- 
arators 


1  pps 


Time 
Inten/al 
Counter 


T 

Ref 


Data 


Figure  34.  Block  diagram  of  a  time  interval  counter  measuring  system. 

This  measurement  method  is  made  practical  by  modem  high-resolution  interpolating  time  interval  counters  that  offer 
10  digit/s  resolution.  The  resolution  is  not  affected  by  the  division  ratio,  which  sets  the  minimum  measurement  time, 
and  determines  how  long  data  can  be  taken  before  a  phase  spillover  occurs  (which  can  be  hard  to  remove  from  a  data 
set).  A  source  having  a  frequency  offset  of  1  x  10"^  can,  for  example,  be  measured  for  only  about  5.8  days  before  a  1 
pps  phase  spillover  occurs  after  being  initially  set  at  the  center.  Drift  in  the  trigger  point  of  the  counter  can  be  a 
limitation  to  this  measurement  method. 

9.2.   Heterodyne  Method 


The  heterodyne  method  mixes  (subtracts)  the  two  sources  being  compared,  and  measures  the  period  of  the  resulting 
audio-frequency  beat  note.  The  measurement  resolution  is  increased  by  the  heterodyne  factor  (the  ratio  of  the  carrier 
to  the  beat  frequency). 
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Figure  35.  Block  diagram  of  a  heterodyne  measuring  system. 

This  heterodyne  technique  is  a  classic  way  to  obtain  high  resolution  with  an  ordinary  period  counter.  It  is  based  on  the 
principle  that  phase  information  is  preserved  in  a  mixing  process.  For  example,  mixing  a  10  MHz  source  against  a 
9.9999  MHz  Hz  offset  reference  will  produce  a  100  Hz  beat  signal  whose  period  variations  are  enhanced  by  a  factor 
of  10  MHz/100  Hz  =  10^  Thus  a  period  counter  with  100  ns  resolution  (10  MHz  clock)  can  resolve  clock  phase 
changes  of  1  ps.  A  disadvantage  of  this  approach  is  that  a  stable  offset  reference  is  required  at  exactly  the  right 
frequency.  Even  worse,  it  can  measure  only  frequency,  requires  a  priori  knowledge  of  the  sense  of  the  frequency 
difference,  and  often  has  dead  time  between  measurements. 

9.3.   Dual  Mixer  Time  Difference  Method 

The  third  method,  in  effect,  combines  the  best  features  of  the  tlrst  two,  using  a  time  interval  counter  to  measure  the 
relative  phase  of  the  beat  signals  fi-om  a  pair  of  mixers  driven  from  a  common  offset  reference: 


^ref 

Figure  36.  Block  diagram  of  a  dual  mixer  time  difference  measuring  system. 

This  dual  mixer  time  difference  (DMTD)  setup  is  arguably  the  most  precise  way  of  measuring  an  ensemble  of  clocks 
all  having  the  same  nominal  frequency.  When  expanded  to  multiple  channels  by  adding  additional  buffer  amplifiers 
and  mixers,  and  time  tagging  the  zero-crossings  of  the  beat  notes  for  each  channel,  this  arrangement  allows  any  two  of 
the  clocks  to  be  intercompared.  The  offset  reference  need  not  be  coherent,  nor  must  it  have  particularly  low  noise  or 
high  accuracy,  because  its  effect  cancels  out  in  the  overall  measurement  process.  For  best  cancellation,  the  zero- 
crossings  should  be  coincident  or  interpolated  to  a  common  epoch.  Additional  counters  can  be  used  to  count  the 
whole  beat  note  cycles  to  eliminate  their  ambiguity,  or  the  zero-crossings  can  simply  be  time  tagged.  The  measuring 
system  resolution  is  determined  by  the  time  interval  counter  or  time-tagging  hardware,  and  the  mixer  heterodyne 
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factor.  For  example,  if  two  5  MHz  sources  are  mixed  against  a  common  5  MHz  to  10  Hz  offset  oscillator  (providing  a 
5  X  10''/10  =  5  X  10'  heterodyne  factor),  and  the  beat  note  is  time  tagged  with  a  resolution  of  100  ns(10  MHz  clock), 
the  measuring  overall  system  resolution  is  10"'' /5  x  10'  =  0.2  ps. 

Multichannel  DMTD  clock  measuring  systems  have  been  utilized  by  leading  national  and  commercial  metrology 
laboratories  for  a  number  of  years  [1-5].  An  early  commercial  version  is  described  in  Reference  [3],  and  a  newer 
technique  is  described  in  Reference  [8].  A  direct  digital  synthesizer  (DDS)  can  be  used  as  the  offset  reference  to 
allow  measurements  to  be  made  at  any  nominal  frequency  within  its  range.  Cross-correlation  methods  can  be  used  to 
reduce  the  DDS  noise.  Instruments  using  those  techniques  are  available  that  automatically  make  both  time  and 
frequency  domain  measurements. 

9.4.  Measurement  Problems  and  Pitfalls 

It  can  be  difficult  to  distinguish  between  a  bad  unit  under  test  and  a  bad  measurement.  When  problems  occur  in  time- 
domain  frequency  stability  measurements,  they  usually  cause  results  that  are  worse  than  expected.  It  is  nearly 
impossible  for  a  measurement  problem  to  give  better  than  correct  results,  and  there  is  considerable  justification  in 
saying  that  the  best  results  are  the  correct  ones.  Two  possible  exceptions  to  this  are  (1)  misinterpretation  of  the  scale 
factor,  and  (2)  inadvertent  coherency  (e.g.,  injection  locking  of  one  source  to  another  due  to  inadequate  isolation. 
Lack  of  stationarity  (changes  in  the  source  itself),  while  not  a  measurement  problem  per  se,  must  also  be  considered. 
In  general,  the  more  devices  available  and  the  more  measurements  being  made,  the  easier  it  is  to  sort  things  out. 

One  common  problem  is  hum  that  contaminates  the  measurements  due  to  ground  loops.  Because  the  measurement 
interval  is  usually  much  longer  than  the  period  of  the  power  line  frequency,  and  not  necessarily  coherent  with  it, 
aliased  "beats"  occur  in  the  frequency  record.  Inspection  of  the  raw  data  can  show  this,  and  the  best  cure  is  often 
isolation  transformers  in  the  signal  leads.  In  fact,  this  is  a  wise  precaution  to  take  in  all  cases. 

All  sorts  of  other  mechanisms  (electrical,  electromagnetic,  magnetic,  thermal,  barometric,  vibrational,  acoustic,  etc.) 
exist  that  can  interfere  with  time  domain  frequency  measurements.  Think  about  all  the  environmental  sensitivities  of 
the  unit  under  test,  and  control  as  many  of  them  as  possible.  Be  alert  to  day-night  and  weekly  cycles  that  indicate 
human  interference.  Stories  abound  about  correlations  between  elevators  moving  and  cars  coming  and  going  ("auto- 
correlations") that  have  affected  clock  measurements.  Think  about  what  can  have  changed  in  the  overall  test  setup. 
Slow  periodic  fluctuations  will  show  up  more  distinctly  in  an  all  tau  (rather  than  an  octave  tau)  stability  plot. 

In  high-precision  measurements,  where  picoseconds  matter  (e.g.  1  x  10""  =  1  ps/1000  seconds),  it  is  important  to 
consider  the  mechanical  rigidity  of  the  test  setup  (e.g.  1  ps  =  0.3  mm).  This  includes  the  electrical  length  (phase 
stability)  of  the  connecting  cables.  Teflon  dielectric  is  an  especially  bad  choice  for  temperature  stability,  while 
foamed  polyethylene  is  much  better.  Even  a  few  degrees  of  temperature  variation  will  cause  the  phase  of  a  high- 
stability  source  to  "breathe"  as  it  passes  through  100  ft  of  coaxial  cable. 

Phase  jumps  are  a  problem  that  should  never  be  ignored.  Examination  of  the  raw  phase  record  is  critical  because  a 
phase  jump  (frequency  impulse)  appears  strangely  in  a  frequency  record  as  a  white  FM  noise  characteristic  [10]. 
Some  large  phase  jumps  are  related  to  the  carrier  period  (e.g.,  a  malfunctioning  digital  frequency  divider). 

It  is  difficult  to  maintain  the  integrity  of  a  measuring  system  over  a  long  period,  but,  as  long  as  the  operating 
conditions  of  the  unit  under  test  and  the  reference  are  undisturbed,  gaps  in  the  data  record  may  be  acceptable.  An 
uninterruptible  power  system  is  indispensable  to  maintain  the  continuity  of  a  long  run. 

9.5.  Measuring  System  Summary 

A  comparison  of  the  relative  advantages  and  disadvantages  of  these  methods  is  shown  in  the  following  table: 
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Table  21.  Comparison  of  time  and  frequency  measurement  methods. 


Method 

Advantages 

Disadvantages 

Divider  &  time  interval 
counter 

Provides  phase  data 

IVlULiC^L  ICbUlUlUJIl 

Covers  wide  range  of  carrier 
frequencies 

Easily  expandable  at  low  cost 

Not  suitable  for  short 
tau 

Mixer  and  period 
counter 

Resolution  enhanced  by 
heterodyne  factor 

No  nhase  data 

Provides  direct  frequency  data 

No  frequency  sense 

Usable  for  short  tau 

Requires  offset 
reference 

Expandable  at  reasonable  cost 

Single  carrier 
frequency 

Dual  mixer  time 
difference 

High  resolution,  low  noise 

Single  carrier 
frequency 

Provides  phase  data 

Offset  reference  noise  & 
inaccuracy  cancels 

Relatively  complex 

No  fixed  reference  channel 

It  is  preferable  to  make  continuous  zero-dead-time  phase  measurements  at  regular  intervals,  and  a  system  using  a 
dual-mixer  time  interval  measurement  is  recommended.  An  automated  high-resolution  multi-channel  clock  (phase) 
measuring  system  with  a  high-performance  (e.g.,  hydrogen  maser)  reference  is  a  major  investment,  but  one  that  can 
pay  off  in  better  productivity.  It  is  desirable  that  the  measurement  control,  data  storage,  and  analysis  functions  be 
separated  to  provide  robustness  and  networked  access  to  the  data.  A  low-noise  reference  not  only  supports  better 
measurement  precision  but  also  allows  measurements  to  be  made  faster  (with  less  averaging). 

9.6.  Data  Format 

A  one-column  vector  is  all  that  is  required  for  a  phase  or  frequency  data  array.  Because  the  data  points  are  equally 
spaced,  no  time  tags  are  necessary.  Nevertheless,  the  use  of  time  tags  is  recommended  (see  section  9.4  below), 
particularly  to  identify  anomalies  or  to  compare  several  sources.  Time  tagging  is  generally  required  for  archival 
storage  of  clock  measurements,  but  a  single  vector  of  extracted  gap-filled  data  is  sufficient  for  analysis.  The 
recommended  unit  for  phase  data  is  seconds,  while  frequency  data  should  be  in  the  form  of  dimensionless  fractional 
frequency.  Double-precision  exponential  ASCII  numeric  format  is  recommended  for  ease  of  reading  into  most 
analysis  software,  with  comma  or  space-delimited  fields  and  one  data  point  per  line.  The  inclusion  of  comments  and 
headers  can  pose  problems,  but  most  software  will  reject  lines  that  start  with  a  or  some  other  non-numeric 
character. 

9.7.  Data  Quantization 

The  phase  or  frequency  data  must  be  gathered  with  sufficient  resolution  to  show  the  variations  of  interest,  and  it  must 
be  represented  with  sufficient  precision  to  convey  those  variations  after  removal  of  fixed  offsets  (see  section  10.1 
below).  Nevertheless,  highly  quantized  data  can  still  contain  useful  information,  especially  after  they  are  combined 
into  longer  averaging  times.  An  example  of  highly  quantized  frequency  data  is  the  random  telegraph  signal  shown 
below.  Although  these  data  have  a  non-Gaussian  amplitude  distribution  (their  histogram  consists  of  two  spikes),  the 
random  occurrences  of  the  two  levels  produce  a  white  FM  noise  characteristic. 
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FREQUENCY  STABILITY 

W  FM  Noise  of  Random  Telegraph  Signal 


Siqno 
tOOe+00 
7.O0e-01 
4,98e-01 
3. 52  e- 01 
2.53e-0i 
'i.85e-0i 
l.32e-01 
9.8Qe-02 
256  6.16e-02 
512  4.6Ce-02 
3.66e-02 
2.5Ge-02 


Freq  Dqla  =  +/— 1  ot  rqndc 


10^  10'  10' 

Aveiaging  Time,  t.  Seconds 


10^ 


Figure  37.  Random  telegraph  signal  as  an  example  of  highly  quantized  frequency  data. 

9.8.  Time  tags 

Time  tags  are  often  associated  witii  phase  or  frequency  data,  and  can  be  usefully  applied  to  the  analysis  of  these  data. 

Time  tags  are  highly  desirable  for  frequency  stability  measurements,  particularly  for  identifying  the  exact  time  of 
some  anomaly.  The  preferred  time  tag  is  the  Modified  Julian  Date  (MJD)  expressed  as  a  decimal  fraction  and 
referenced  to  UTC.  Based  on  the  astronomical  Julian  Date,  the  number  of  days  since  noon  on  January  1,  4713  BC, 
the  MJD  is  the  Julian  Date  -  2  4000  000.5.  It  is  widely  used,  purely  numeric,  can  have  any  required  resolution,  is 
easily  converted  to  other  formats,  is  non-ambiguous  over  a  two-century  (1900  to  2099)  range,  and  is  free  from 
seasonal  (daylight  saving  time)  discontinuities 

Analysis  software  can  easily  convert  the  MJD  into  other  formats  such  as  year,  month,  day,  hour,  minute,  and  second. 
The  MJD  (including  the  fractional  day)  can  be  obtained  from  the  C  language  timeo  function  by  dividing  its  return 
value  by  86  400  and  adding  40  587. 

9.9.  Archiving  and  Access 


There  is  no  standard  way  to  archive  and  access  clock  data.  For  some  purposes,  it  is  sufficient  to  simply  save  the  raw 
phase  or  frequency  data  to  a  file,  identifying  it  only  by  the  file  name.  At  the  other  extreme,  multichannel  clock 
measuring  systems  may  require  an  elaborate  database  to  store  a  large  collection  of  data,  keep  track  of  the  clock 
identities  and  transactions,  provide  security  and  robust  data  integrity,  and  serve  the  archived  data  via  a  network.  It 
may  also  be  necessary  to  integrate  the  clock  data  with  other  information  (e.g.,  temperature)  from  a  data  acquisition 
system. 
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10     Analysis  Procedure 


A  frequency  stability  analysis  can  proceed  along  several  paths,  as  the  circumstances  dictate.  Nevertheless,  the 
following  example  shows  a  typical  analysis  flow.  Using  simulated  data  for  a  high-stability  rubidium  frequency 
standard,  the  purpose  of  the  analysis  is  to  characterize  the  noise  in  the  presence  of  an  outlier,  large  frequency  offset 
and  significant  drift. 


Table  22.  An  example  of  a  typical  analysis  flow. 


Step 


Description 


Plot 


Open  and  examine  a 
phase  data  file.  The 
phase  data  is  just  a 
ramp  with  slope 
corresponding  to 
frequency  offset. 


0000 
9998 
7998 
6998 
5998 
4998 
3998 
2998 
1998 
0997 
9997 
8997 


000000 
730257 
905261 
692150 
738512 
876276 
361914 
336122 
367234 
856792 
748960 
326982 
etc . 


OOOOOe+00 
41449e-07 
85009e-06 
98003e-06 
09537e-06 
63997e-06 
40859e-06 
16789e-06 
54638e-06 
57264e-06 
24524e-06 
42008e-06 


PHASE  DATA 

Simulated  RAF3  Phase  Data 


10.0    12J    15.0    17.5    20.0    22.5    25.0    27.5  30.0 

Time,  Days 


Convert  the  phase  data 
to  frequency  data  and 
examine  it.  An 
obvious  outlier  exists 
that  must  be  removed 
to  continue  the 
analysis. 

Visual  inspection  of 
data  is  an  important 
preprocessing  step! 

Analyst  judgment  may 
be  needed  for  less 
obvious  outliers. 


FREQUENCY  DATA 

Simulated  RAFS  Frequency  Data 


2.5     5.0  7.5 


10.0    12.5    15.0    17.5    20.0    22.5    25.0    27J  30.0 

Time,  Days 
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In  an  actual  analysis, 
one  should  try  to 
determine  the  cause  of 
the  outlier.  The 
frequency  spike  of 
IxlO  ''  corresponds  to  a 
phase  step  of  900  ns 
over  a  single  900-s 
measurement  interval, 
nine  10  MHz  carrier 
cycles.  Data  taken  at  a 
higher  rate  would  help 
to  determine  whether 
the  anomaly  happened 
instantaneously  or  over 
some  finite  period. 
Timetags  can  help  to 
relate  the  outlier  to 
other  external  events. 


PHASE  DATA 

Phase  Slep  at  Frequency  Outlier 


-«0 
-500 
-600 
.700 
800 
-900 
lOOO 


435.0  1436.0  1437.0  1438.0  1439.0  1440.0  1441.0  1442.0  1443.0  1444.0  1445.0 

Data  Point 


Remove  the  outlier. 
The  noise  and  drift  are 
now  visible.  A  line 
shows  a  linear  fit  to  the 
frequency  data,  which 
appears  to  be  quite 
appropriate. 


FREQUENCY  DATA 

Frequency  Outlier  Removed 


2J     5.0      7.5     10.0    12.5    15.0    17J    20.0    22.5    25.0    27.5    30  0 

Time,  Days 


Remove  the  frequency 
offset  from  the  phase 
data.  The  resulting 
quadratic  shape  is  due 
to  the  frequency  drift. 
One  can  just  begin  to 
see  phase  fluctuations 
around  the  quadratic  fit 
to  the  phase  data. 


PHASE  DATA 

Frequency  Oltset  Removed 


2.5      5.0      7J     10.0     12J     15.0     17.5    20.0    22.5    25.0    27.5  30.0 

Time,  Days 
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Remove  the  frequency 
drift,  leaving  the  phase 
residuals  for  noise 
analysis,  which  is  now 
clearly  visible. 

Some  experience  is 
needed  to  interpret 
phase  data  like  these. 
Remember  that 
frequency  corresponds 
to  the  slope  of  the 
phase,  so  the  frequency 
is  lowest  near  the  end 
of  the  record,  where  the 
phase  slope  is  the  most 
negative.  


PHASE  DATA 

Frequency  Drift  Removed 


7.5     10.0    12.3    15.0    17J    20.0    22.5  25.0 

Time,  Days 


Convert  the  phase 
residuals  to  frequency 
residuals. 

Alternatively,  remove 
the  frequency  drift 
from  the  frequency 
data  of  Step  #3.  There 
are  subtle  differences 
in  removing  the  linear 
frequency  drift  as  a 
quadratic  fit  to  the 
phase  data  compared 
with  removing  it  as  a 
linear  fit  to  the 
frequency  data 
(different  noise  models 
apply).  Other  drift 
models  may  be  more 
appropriate.  Analyst 
judgment  is  needed  to 
make  the  best  choices. 


FREQUENCY  DATA 

Frequency  Onset  and  Drift  Removed 


10.0    12.5    15.0    17J    20.0    22.5  25.0 

Time,  Days 


Perform  a  stability 
analysis  using  the 
overlapping  Allan 
deviation.  The  results 
show  white  FM  noise 
at  short  averaging 
times  (i"''^  slope)  ,  and 
flicker  FM  noise  at 
longer  tau  (t°  slope), 
both  at  the  simulated 
levels  shown  in  the 
annotations  of  the  first 
plot.  


FREQUENCY  STABILITY 

 Frequency  Drift  Removed  


2.<7e- 
1.95e- 
1  55e- 


2-30e  +  05  1-18 


lO'  10*  lO' 

Averaging  Time,  T,^  Seconds 
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10.1.  Data  Precision 


There  are  relatively  few  numerical  precision  issues  relating  to  the  analysis  of  frequency  stabilit>'  data.  One  exception, 
however,  is  phase  data  for  a  highly  stable  frequency  source  having  a  relatively  large  frequency  offset.  The  raw  phase 
data  will  be  essentialh  a  straight  line  (representing  the  frequency  offset),  and  the  instability  information  is  contained 
in  the  small  deviations  from  the  line.  A  large  number  of  digits  must  be  used  unless  the  frequency  offset  is  removed  by 
subtracting  a  linear  term  from  the  raw  phase  data.  Similar  considerations  apply  to  the  quadratic  phase  term  (linear 
frequency  drift).  Man\  frequency  stability,  measures  involve  averages  of  tlrst  or  second  differences.  Thus,  while  their 
numerical  precision  obviously  depends  upon  the  variable  digits  of  the  data  set,  there  is  little  error  propagation  in 
forming  the  summan.  statistics. 


10.2.  Preprocessing 

Preprocessing  of  the  measurement  data  is  often  necessary  before  the  actual  analysis  is  performed,  which  may  require 
data  averaging,  or  removal  of  outliers,  frequency  offset,  and  drift. 


Phase  data  may  be  converted  to  frequencv  data,  and  vice  versa.  Phase  and  frequency  data  can  be  combined  for  a 
longer  averaging  time.  Frequenc\  offset  ma>  be  removed  from  phase  data  by  subtracting  a  line  determined  b\  the 
average  of  the  first  differences,  or  by  a  least  squares  linear  fit.  An  offset  may  be  removed  from  frequency  data  by 
normalizing  it  to  have  an  average  value  of  zero.  Frequency  drift  may  be  removed  from  phase  data  by  a  least  squares 
or  three-point  quadratic  fit,  or  by  subtracting  the  average  of  the  second  differences.  Frequency  drift  may  be  removed 
from  frequency  data  by  subtracting  a  least-squares  linear  fit.  by  subtracting  a  line  determined  by  the  first  differences 
or  by  calculating  the  drift  from  the  difference  between  the  two  halves  of  the  data.  The  latter,  called  the  bisection  drift, 
is  equivalent  to  the  three-point  fit  for  phase  data.  Other  more  specialized  log  and  diffusion  models  may  also  be  used. 
The  latter  are  particularh  useful  to  describe  the  stabilization  of  a  frequency  source.  In  general,  the  objective  is  to 
remove  as  much  of  the  deterministic  behavior  as  possible,  obtaining  random  residuals  for  subsequent  noise  analysis. 


10.3.  Gaps,  Jumps,  and  Outliers 


It  is  common  to  have  gaps  and  outliers  in  a  set  of  raw  frequency  stability  data. 
Missing  or  erroneous  data  ma\  occur  due  to  power  outages,  equipment 
malfunctions,  and  interference.  For  long-term  tests,  it  may  not  be  possible  or 
practical  to  repeat  the  run.  or  otherwise  avoid  such  bad  data  points.  Usually  the 
reason  for  the  gap  or  outlier  is  known,  and  it  is  particular!)  important  to  e.xplain  all 
phase  discontinuities.  Plotting  the  data  will  often  show  the  bad  points,  which  may 
have  to  be  removed  before  doing  an  analysis  to  obtain  meaningful  results. 


Gaps,  jumps  and  outliers  can 
occur  in  frequency 
measurements  and  they  must 
be  handled  before  performing 
a  stability  analysis.  Methods 
are  available  to  fill  gaps  and 
to  correct  for  outliers  in  a 
consistent  manner. 


Frequency  outliers  are  found  by  comparing  each  data  point  with  the  median  value  of 

the  data  set  plus  or  minus  some  multiple  of  the  absolute  median  deviation.  These  median  statistics  are  more  robust 
because  they  are  insensitive  to  the  size  of  the  outliers.  Outliers  can  be  replaced  by  gaps  or  filled  w  ith  interpolated 
values. 


Frequency  jumps  can  also  be  a  problem  for  stabilitv  analysis.  Their  occurrence  indicates  that  the  statistics  are  not 
stationar>  ,  and  it  may  be  necessar>  to  divide  the  data  into  portions  and  analyze  them  separately. 


Gaps  and  outliers  can  occur  in  clock  data  due  to  problems  with  the  measuring  system  or  the  frequency  source  itself 
Like  death  and  taxes,  gaps  and  outliers  can  be  avoided  but  not  eliminated. 
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10.4.  Gap  Handling 


Gaps  should  be  included  to  maintain  the  proper  implied  time  interval  between  measurements,  and  a  value  of  zero  (0) 
is  often  used  to  denote  a  gap.  For  phase  data,  zero  should  be  treated  as  valid  data  if  it  is  the  first  or  last  point.  For 
fractional  frequency  data,  valid  data  having  a  value  of  zero  can  be  replaced  by  some  very  small  value  (e.g.,  le-99). 
Many  analysis  functions  can  produce  meaningful  results  for  data  having  gaps  by  simply  skipping  those  points  that 
involve  a  gap.  For  example,  in  the  calculation  of  the  Allan  variance  for  frequency  data,  if  either  of  the  two  points 
involved  in  the  first  difference  is  a  gap,  that  Allan  variance  pair  is  skipped  in  the  summation. 

Gaps  may  be  filled  in  phase  or  frequency  data  by  replacing  them  with  interpolated  values,  by  first  removing  any 
leading  and  trailing  gaps,  and  then  using  the  two  values  immediately  before  and  after  any  interior  gaps  to  determine 
linearly  interpolated  values  within  the  gap. 

A  zero  value  in  fractional  frequency  data  can  also  occur  as  the  result  of  the  conversion  of  two  equal  adjacent  phase 
data  points  (perhaps  because  of  limited  measurement  resolution),  and  the  value  should  be  adjusted  to,  say,  le-99  to 
distinguish  it  from  a  gap. 

10.5.  Uneven  Spacing 

Unevenly  spaced  phase  data  can  be  handled  if  they  have  associated  time  tags  by  using  the  individual  time  tag  spacing 
when  converting  it  to  frequency  data.  Then,  if  the  tau  differences  are  reasonably  small,  the  data  may  be  analyzed  by 
using  the  average  time  tag  spacing  as  the  analysis  tau,  in  effect  placing  the  frequency  data  on  an  average  uniform  grid. 
While  completely  random  data  spacing  is  not  amenable  to  this  process,  tau  variations  of  ±10  %  will  yield  reasonable 
results  as  long  as  the  exact  intervals  are  used  for  the  phase  to  frequency  conversion. 

10.6.  Analysis  of  Data  with  Gaps 

Care  must  be  taken  when  analyzing  the  stability  of  data  with  missing  points  and/or  gaps.  Missing  points  can  be  found 
by  examining  the  time  tags  associated  with  the  data,  and  gaps  can  then  be  inserted  as  placeholders  to  maintain  equally 
spaced  data.  Similarly,  outliers  can  be  replaced  with  gaps  for  the  same  reason.  These  gaps  can  span  multiple  points. 
Some  analysis  processes  can  be  performed  with  data  having  gaps  by  skipping  over  them,  perhaps  at  some  speed 
penalty,  but  other  calculations  cannot  be.  It  is  therefore  often  necessary  to  replace  the  gaps  with  interpolated  values. 
Those  points  are  not  real  data,  however,  and,  if  there  are  many  of  them,  the  results  will  be  suspect.  In  these  cases, 
judgment  is  needed  to  assure  a  credible  result.  It  may  be  more  prudent  to  simply  analyze  a  gap-free  portion  of  the 
data. 

10.7.  Phase-Frequency  Conversions 

Phase  to  frequency  conversion  is  straightforward  for  data  having  gaps.  Because  two  phase  points  are  needed  to 
determine  each  frequency  point  (as  the  difference  between  the  phase  values  divided  by  their  tau),  a  single  phase  gap 
will  cause  two  frequency  gaps,  and  a  gap  of  N  phase  points  causes  N  +  1  frequency  gaps. 

Conversion  from  frequency  to  phase  is  more  problematic  because  of  the  need  to  integrate  the  frequency  data.  The 
average  frequency  value  is  used  to  calculate  the  phase  during  the  gap,  which  can  cause  a  discontinuity  in  the  phase 
record.  Analysis  of  phase  data  resulting  from  the  conversion  of  frequency  data  having  a  large  gap  is  not 
recommended.  ♦ 
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10.8.  Drift  Analysis 


Drift  analysis  fianctions  generally  perform  well  for  data  having  gaps,  provided  that  missing  data  are  represented  by 
gaps  to  maintain  a  regular  time  sequence. 

10.9.  Variance  Analysis 

Variance  analysis  functions  can  include  provisions  for  handling  gaps.  Some  of  these  functions  yield  satisfactory 
results  in  all  cases,  while  others  have  speed  limitations,  or  provide  unsatisfactory  results  for  data  having  large  gaps. 
The  latter  is  most  apparent  at  longer  averaging  times  where  the  averaging  factor  is  comparable  to  the  size  of  the  gap. 
The  speed  limitations  are  caused  by  more  complex  gap  checking  and  frequency  averaging  algorithms,  while  the  poor 
results  are  associated  with  the  total  variances  for  which  conversion  to  phase  data  is  required.  In  all  cases,  the  results 
will  depend  on  coding  details  included  in  addition  to  the  basic  variance  algorithm.  Filling  gaps  can  often  help  for  the 
total  variances.  Two  general  rules  apply  for  the  variance  analysis  of  data  having  large  gaps:  (1)  use  unconverted 
phase  data,  and  (2)  check  the  results  against  the  normal  Allan  deviation  (which  has  the  simplest,  fastest  gap  handling 
ability). 

10.10.  Spectral  Analysis 

Gap  filling  in  spectral  analysis  ftjnctions  can  affect  the  low  frequency  portion  of  the  spectrum. 

10.11.  Outlier  Recognition 

The  median  absolute  deviation  (MAD)  is  recommended  as  its  means  of  outlier  recognition.  The  MAD  is  a  robust 
statistic  based  on  the  median  of  the  data.  It  is  the  median  of  the  scaled  absolute  deviations  of  the  data  points  from 
their  median  value,  defined  as  MAD  =  Median  {  |  y(i)  -  m  |  /  0.6745  },  where  m  =  Median  {  y(i)  },  and  the  factor 
0.6745  makes  the  MAD  equal  to  the  standard  deviation  for  normally  distributed  data.  Each  frequency  data  point,  y(i), 
is  compared  with  the  median  value  of  the  data  set,  m,  plus  or  minus  the  desired  multiple  of  the  MAD. 

While  the  definition  of  an  outlier  is  somewhat  a  matter  of  judgment,  it  is  important  to  find  and  remove  such  points  in 
order  to  use  the  rest  of  the  data,  based  on  their  deviation  from  the  median  of  the  data,  using  a  deviation  limit  in  terms 
of  the  median  absolute  deviation  (a  5-sigma  limit  is  common).  This  is  a  robust  way  to  determine  an  outlier,  which  is 
then  replaced  by  a  gap.  An  automatic  outlier  removal  algorithm  can  iteratively  apply  this  method  to  remove  all 
outliers,  which  should  be  an  adjunct  to,  and  not  a  substitute  for,  visual  inspection  of  the  data. 

It  is  important  to  explain  all  outliers,  thereby  determining  whether  they  are  due  to  the  measurement  process  or  the 
device  under  test.  An  important  first  step  is  to  correlate  the  bad  point  with  any  external  events  (e.g..  power  outages, 
equipment  failures,  etc.)  that  could  account  for  the  problem.  Failures  of  the  measurement  system,  frequency 
reference,  or  environmental  control  are  often  easier  to  identify  if  multiple  devices  are  under  test.  Obviously,  a 
common  gap  in  all  measurement  channels  points  to  a  failure  of  the  measurement  system,  while  a  common  change  in 
all  measurement  readings  points  to  a  reference  problem.  Auxiliary  information  such  as  monitor  data  can  be  a  big  help 
in  determining  the  cause  of  outliers.  A  log  of  all  measurement  system  events  should  be  kept  to  facilitate  outlier 
identification.  Discrete  phase  jumps  are  a  particular  concern,  and,  if  they  are  related  to  the  RF  carrier  frequency,  may 
indicate  a  missing  cycle  or  a  problem  with  a  digital  divider.  A  phase  jump  will  correspond  to  a  frequency  spike  with  a 
magnitude  equal  to  the  phase  change  divided  by  the  measurement  interval.  Such  a  frequency  spike  will  produce  a 
stability  record  that  appears  to  have  a  (large  magnitude)  white  FM  noise  characteristic,  which  can  be  a  source  of 
confusion. 
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10.12.  Data  Plotting 

Data  plotting  is  often  the  most  important  step  in  the  analysis  of  frequency  stability.  Visual  inspection  can  provide  vital 
insight  into  the  results,  and  is  an  important  preprocessor  before  numerical  analysis.  A  plot  also  shows  much  about  the 
validity  of  a  curve  fit. 

Phase  data  are  plotted  as  line  segments  connecting  the  data  points.  This  presentation  properly  conveys  the  integral 
nature  of  the  phase  data.  Frequency  data  are  plotted  as  horizontal  lines  between  the  frequency  data  points.  This  shows 
the  averaging  time  associated  with  the  frequency  measurement,  and  mimics  the  analog  record  from  a  frequency 
counter.  As  the  density  of  the  data  points  increases,  there  is  essentially  no  difference  between  the  two  plotting 
methods.  Missing  data  points  are  shown  as  gaps  without  lines  connecting  the  adjacent  points. 

10.13.  Variance  Selection 

It  is  the  user's  responsibility  to  select  an  appropriate  variance  for  the  stability  analysis.  The  overlapping  Allan 
variance  is  recommended  in  most  cases,  especially  where  the  frequency  drift  is  small  or  has  been  removed.  The  Allan 
variance  is  the  most  widely  used  time-domain  stability  measure,  and  the  overlapping  form  provides  better  confidence 
than  the  original  "normal"  version.  The  total  and  TheoH  variance  can  be  used  for  even  better  confidence  at  large 
averaging  factors  (but  at  the  expense  of  longer  computation  time).  The  modified  Allan  variance  is  recommended  to 
distinguish  between  white  and  flicker  PM  noise,  and,  again,  a  total  form  of  it  is  available  for  better  confidence  at  long 
tau.  The  time  variance  provides  a  good  measure  of  the  time  dispersion  of  a  clock  due  to  noise,  while  MTIE  measures 
the  peak  time  deviations.  TIE  rms  can  also  be  used  to  assess  clock  performance,  but  TVAR  is  generally  preferred. 
Finally,  the  overlapping  Hadamard  variance  is  recommended  over  its  normal  form  for  analyzing  stability  in  the 
presence  of  divergent  noise  or  frequency  drift.  In  all  cases,  the  results  are  reported  in  terms  of  the  deviations. 

The  choice  of  tau  interval  depends  mainly  on  whether  interference  mechanisms  are  suspected  that  cause  the  stability 
to  vary  periodically.  Normally,  octave  or  decade  spacing  is  used  (the  former  has  even  spacing  on  a  log-log  plot,  while 
the  latter  provides  tau  multiples  of  ten).  The  all  tau  option  can  be  useful  as  a  form  of  spectral  analysis  to  detect  cyclic 
disturbances  (such  as  temperature  cycling). 

10.14.  Three-Cornered  Hat 


Any  frequency  stability  measurement  includes  noise  contributions  from  both  the 
device  under  test  and  the  reference.  Ideally,  the  reference  noise  would  be  low 
enough  that  its  contribution  to  the  measurement  is  negligible.  Or,  if  the  noise  of  the 
reference  is  known,  it  can  be  removed  by  subtracting  its  variance.  A  special  case  is 
that  of  two  identical  units  where  half  of  the  measured  variance  comes  from  each, 
and  the  measured  deviation  can  be  corrected  for  one  unit  by  dividing  it  by  V2. 
Otherwise,  it  may  be  useful  to  employ  the  so-called  "three-cornered  haf  method  for 
determining  the  variance  of  an  individual  source.  Given  a  set  of  three  pairs  of 
measurements  for  three  independent  frequency  sources  a,  b,  and  c  whose  variances 


It  is  sometimes  necessary  to 
determine  the  individual 
noise  contributions  of  the  two 
sources  that  contribute  to  a 
variance  measurement. 
Methods  exist  for  doing  so  by 
using  the  results  of  multiple 
measurements  to  estimate  the 
variance  of  each  source. 
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The  individual  variances  may  be  determined  by  the  expressions 
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Although  useful  for  determining  the  individual  stabilities  of  units  having  similar  performance,  the  method  may  fail  by 
producing  negative  variances  for  units  that  have  widely  differing  stabilities,  if  the  units  are  correlated,  or  for  which 
there  are  insufficient  data.  The  three  sets  of  stability  data  should  be  measured  simultaneously.  The  three-cornered  hat 
method  should  be  used  with  discretion,  and  it  is  not  a  substitute  for  a  low  noise  reference.  It  is  best  used  for  units 
having  similar  stability  (e.g.,  to  determine  which  unit  is  best).  Negative  variances  are  a  sign  that  the  method  is  failing 
(because  it  was  based  on  insufficient  measurement  data,  or  because  the  units  under  test  have  disparate  or  correlated 
stability ).  This  problem  is  most  likely  to  arise  at  long  tau. 

The  three-cornered  hat  function  may  be  used  to  correct  a  stability  measurement  for  the  noise  contribution  of  the 
reference,  as  shown  in  Figure  38. 


Figure  38.  Three-cornered  hat  function. 

The  Unit  Under  Test  (UUT),  denoted  as  source  A,  is  measured  against  the  reference,  denoted  by  B  and  C,  by  identical 
stability  data  files  A-B  and  A-C.  The  reference  is  measured  against  itself  by  stability  data  file  B-C,  which  contains  the 
a  priori  reference  stability  values  multiplied  by  V2. 

An  example  of  the  use  of  the  three-cornered  hat  function  to  correct  stability  data  for  reference  noise  is  shown  below. 
Simulated  overlapping  Allan  deviation  stability  data  for  the  unit  under  test  versus  the  reference  were  created  by 
generating  and  analyzing  512  points  of  frequency  data  with  tau  =  1  s  and  ay(l  s)  =  le-1 1.  The  resulting  stability  data 
are  shown  in  the  following  table. 
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Table  23.  Stability  data  for  unit  under  test  versus  reference. 


Tau 

# 

Sigma 

Min  Sigma 

Max  Sigma 

1  .OOOe+00 

51 1 

9.448e-12 

9.108e-12 

9.830e-l  2 

2. OOOe+00 

509 

7.203e-12 

6.923e-12 

7.520e-l  2 

4. OOOe+00 

50:? 

5.075e-12 

4.826e-12 

5.367e-12 

8. OOOe+00 

497 

3.27:)e-l  2 

j>.058e-12 

3.546e-12 

1.600e+01 

481 

2.370e-12 

2.157e-12 

2.663e-12 

3.200e+01 

449 

1.854e-12 

1.720e-12 

2.025e-12 

6.400e+01 

385 

l.269e-12 

1.147e-12 

1.441e-12 

1.280e+02 

257 

5.625e-13 

4.820e-13 

7.039e-13 

A  similar  stability  file  is  used  for  the  reference.  Since  it  represents  a  measurement  of  the  reference  against  itself,  the 
Allan  deviations  of  the  reference  source  are  multiplied  by  V2.  Simulated  overlapping  Allan  deviation  stability  data  for 
the  reference  versus  the  reference  was  created  by  generating  512  points  of  frequency  data  with  tau  =  1  second  and 
ay(l)=  1.414e-12. 


Table  24.  Stability  data  for  unit  under  test  versus  reference. 


Tau 

# 

Sigma 

Min  Sigma 

Max  Sigma 

1  .OOOe+00 

511 

1.490e-12 

1.436e-12 

1.550e-12 

2.000e+00 

509 

1.025e-12 

9.854e-13 

1.070e-12 

4.000e+00 

505 

7.631e-13 

7.257e-13 

8.070e-13 

8.000e+00 

497 

5.846e-13 

5.458e-13 

6.329e-13 

1.600e+01 

481 

3.681e-13 

3.349e-13 

4.135e-13 

3.200e+01 

449 

2.451e-13 

2.152e-13 

2.924e-13 

6.400e+01 

385 

1.637e-13 

1.368e-13 

2.173e-13 

1.280e+02 

257 

1.360e-13 

1.058e-13 

2.285e-13 
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Figure  39.  Corrected  UUT  and  reference  stabilities. 

Here  the  reference  stability  is  about  1  x  10"'"r'  ",  and  the  corrected  UUT  instability  is  slightly  less  than  the 
uncorrected  values.  Note  that  the  B  and  C  columns  of  corrected  stability  values  both  represent  the  reference  source. 
Appropriate  use  of  the  three-cornered  hat  method  to  correct  stability  measurements  for  reference  noise  applies  where 
the  reference  stability  is  between  three  to  10  times  better  than  that  of  the  unit  under  test.  The  correction  is  negligible 
for  more  the  latter  (see  above),  and  has  questionable  confidence  for  less  than  the  former  (and  a  better  reference  should 
be  used). 

The  error  bars  of  the  individual  variances  may  be  set  using  statistics  by  first  determining  the  reduced  number  of 
degrees  of  freedom  associated  with  the  three-cornered  hat  process  [11,  12].  The  fraction  of  remaining  degrees  of 
freedom  for  unit  i  as  a  result  of  performing  a  three-cornered  hat  instead  of  measuring  against  a  perfect  reference  is 
given  by: 


r  = 


2-cr,  -(-cr„ 


1  7  7 


(70) 


The  ratio  of  the  number  of  degrees  of  freedom  is  0.4  for  three  units  having  the  same  stability,  independent  of  the 
averaging  time  and  noise  ty  pe. 

The  three-cornered  hat  technique  can  be  extended  to  M  clocks  (subject  to  the  same  restriction  against  negative 
variances)  by  using  the  expression 
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where 
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and  the  o^,,  are  the  measured  Allan  variances  for  clock  i  versus]  at  averaging  time  t.  Using  o^,,  =  0  and  a2,j=  a^j,,  we  can 
easily  write  closed-form  expressions  for  the  separated  variances  from  measurements  of  M  clocks. 
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10.15.  Reporting 

The  results  of  a  stability  analysis  are  usually  presented  as  a  combination  of  textual,  tabular,  and  graphic  forms.  The 
text  describes  the  device  under  test,  the  test  setup,  and  the  methodology  of  the  data  preprocessing  and  analysis,  and 
summarizes  the  results.  The  assumptions  and  analysis  choices  should  be  documented  so  that  the  results  could  be 
reproduced.  The  report  often  includes  a  table  of  the  stability  statistics.  Graphical  presentation  of  the  data  at  each 
stage  of  the  analysis  is  generally  the  most  important  aspect  of  presenting  the  results.  For  example,  these  are  often  a 
series  of  plots  showing  the  phase  and  frequency  data  with  an  aging  fit,  phase  and  frequency  residuals  with  the  aging 
removed,  and  stability  plots  with  noise  fits  and  error  bars.  Plot  titles,  subtitles,  annotations  and  inserts  can  be  used  to 
clarify  and  emphasize  the  data  presentation.  The  results  of  several  stability  runs  can  be  combined,  possibly  along  with 
specification  limits,  into  a  single  composite  plot.  The  various  elements  can  be  combined  into  a  single  electronic 
document  for  easy  printing  and  transmittal. 
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11     Case  studies 


This  section  contains  several  case  studies  to  further  illustrate  methodologies  of  frequency  stability  analysis. 

11.1.  Flicker  Floor  of  a  Cesium  Frequency  Standard 


The  purpose  of  this  analysis  is  to  examine  the  flicker  floor  of  a  commercial  cesium  beam  tube  frequency  standard.  An 
instrument  of  this  kind  can  be  expected  to  display  a  white  FM  noise  characteristic  out  to  long  averaging  times.  At 
some  point,  however,  the  unit  will  typically  "flicker  out",  with  its  stability  plot  flattening  out  to  a  flicker  FM  noise 
characteristic,  beyond  which  its  stability  will  no  longer  improve.  Determining  this  point  requires  a  lengthy  and 
expensive  run,  and  it  is  worthwhile  to  use  an  analytical  method  that  provides  the  best  information  at  long  averaging 
factors.  The  effort  of  a  more  elaborate  analysis  is  far  easier  than  extending  the  measurement  time  by  weeks  or 
months. 


Table  25.  Flicker  floor  of  a  cesium  frequency  standard. 


This  is  the  five-month 
frequency  record  for  the 
unit  under  test.  It  has 
already  been  "'cleaned 
up"  for  any  missing 
points,  putting  in  gaps  as 
required  to  provide  a 
time-continuous  data  set. 
The  data  look  very 
"white"  (see  Section  6.1), 
the  frequency  offset  is 
very  small  (+7.2  x  10"'^), 
and  there  is  no  apparent 
drift  (+1.4  X  10"'^/day). 
Overall,  this  appears  to 
be  a  very  good  record. 


FREQUENCY  DATA 

Cesium  Beam  Frequency  Standard 


An  overlapping  Allan 
deviation  plot  shows  a 
white  FM  noise  level  of 
about  8.2  X  10'' V"-  out 
to  about  5  days,  where 
the  stability  levels  off  at 
about  1.2  X  10"'^  While 
this  is  very  respectable 
behavior,  one  wonders 
what  the  stability  actually 
is  at  the  longer  averaging 
times.  But  to  gain 
meaningful  confidence 
using  ADEV  there,  the 
run  would  have  to  be 
extended  by  several 
months. 


FREQUENCY  STABILITY 

Cesium  Beam  Frequency  Standard 


10^ 


10^  lO*  lO' 

Averaging^  Time,  x,^  Seconds 
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The  total  deviation  can 
provide  better  confidence 
at  the  longer  averaging 
times.  It  seems  to 
indicate  that  the  stability 
continues  to  improve  past 
10  days,  where  is  drops 
below  1  X  10'"^  but  the 
results  are  not 
conclusive. 


FREQUENCY  STABILITY 

Cesium  Beam  Frequency  Standard 


5i2e*04  3  72e- 

i.02e  +  05  2.5:5e- 

+  05  I59e-K 

)  lOe  +  05  1.24e-1i 

i.I9e-05  1.19e-l4 

.e-06  i.i8e-ii 

3.28e-^06  6  95e- 


lO^ 


10^  10"  10^  10* 

Averaging  Tinie^  J,^  Seconds 


FREQUENCY  STABILITY 

Cesium  Beam  Frequency  Standard 


Theol  can  provide  even 
better  long-term 
information,  at  the 
expense  of  a  longer 
calculation  time.  This 
Theol  plot  seems  to 
show  even  more  clearly 
that  the  stability 
continues  to  improve  at 
longer  averaging  times, 
well  into  the  pplO'^ 
range.  It  assumes  white 
FM  noise  (no  Theol  bias 
removal).  


10^  10''  lO' 

Averaging  Tim&,  x,^  Seconds 


A  TheoH  analysis  with 
automatic  bias  removal 
combines  AVAR  in  short 
and  mid-term  averaging 
times  and  TheoBR  in 
long  term.  TheoH  does 
not  require  any  explicit 
knowledge  of  the  type  of 
noise.  It  appears  that  one 
noise  prevails  in  this  data 
run,  in  which  case  Theol 
(shown  above)  bias  is 
detected  as  1 .676, 
intermediate  between 
white  and  flicker  FM 
noises  for  medium  to 
large  averaging  times. 
Thus,  Theol  results  are 
essentially  the  same  as 
TheoH. 


FREQUENCY  STABILITY 

Cesium  Beam  Frequency  Standard 


ITheol  Bios  =  1  6761  . 

2048 
4095 
8192 


4  69e-14 
J.35e-14 
2.19e-K 
i.53e-K 
i.22e-M 
l.23e-H 


10^ 


10^  lO"  lO' 

Averaging  Time,  t.  Seconds 


We  can  conclude  that  this  cesium  beam  frequency  standard  reaches  a 
stability  slightly  better  than  1  x  lO"''*  at  an  averaging  time  on  the  order  of  1 
month. 
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11.2.  Hadamard  Variance  of  a  Source  with  Drift 


Table  26.  Hadamard  variance  of  a  source  with  drift. 


These  frequency  data 
simulate  a  typical  rubidium 
frequency  standard  (RFS) 
with  a  combination  of  white 
and  flicker  FM  noise,  plus 
significant  frequency  drift. 


FREQUENCY  DATA 

Simulated  RFS  Data 


10.0    12.5    15.0    17.5    20.0   22.5   25.0    27.5  30.0 

Time,  Days 


If  an  Allan  variance  analysis 
is  performed  directly  on 
these  data  without  drift 
removal,  the  stability  at  the 
longer  averaging  times  is 
degraded.  In  particular,  the 
stability  plot  has  a  t"'  slope 
beyond  10^  seconds  that 
corresponds  to  the  drift  (i.e., 
about  1  X  10"'-  at  5  days). 


FREQUENCY  STABILITY 

Slumulaled  RFS  Stablllly 


10>  10'  10* 

Averaging  Time,  T,  Seconds 


If  the  linear  frequency  drift 
is  removed  before 
performing  the  AVAR 
stabilit>  analysis,  the 
stability  plot  shows  a  white 
FM  noise  (t""") 
characteristic  changing  to 
flicker  FM  noise  (t°)  at 
longer  averaging  times.  It  is 
usually  best  to  use  a  stability 
plot  only  to  show  the  noise, 
and  analyzing  and  removing 
the  drift  separately.  


FREQUENCY  STABILITY 

Simulated  RFS  Data 


On 


-cu  Sigma 

9.00e^02  5.50e-lJ 

ISOe^Oi  4,146-13 

3  506^03  3  47e-13 
7  20e-03  3.13e-13 

■ 44e-04  2 75e-13 

2  386^04  2  32e-13 

5  76e-04  2  78e-i3 

:.I5e  +  05  2,42e-:3 

2  30e^05  2.52e-13 

4  61e-05  2.2Ie-I3 


10^  10"  lo' 

Averaging  Tinie,  T,  Seconds 
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The  Hadamard  variance  is 
insensitive  to  linear 
frequency  dritt.  It  can 
therefore  be  used  to  perform 
a  stability  analysis  on  the 
original  data  without  first 
having  to  remove  the  drift. 
The  HDEV  results  are 
essentially  identical  to  those 
of  the  drift-removed  ADEV. 
This  can  be  more 
convenient,  especially  when 
analyzing  a  combination  of 
sources  having  differing 
amounts  of  drift  (e.g., 
cesium  and  rubidium  units). 


FREQUENCY  STABILITY 

simulated  RFS  Data 


9  00e^02  5.43e-l3 

1,806^03  4.00e-i3 

3  60et03  3.216-13 

7  20e*03  2.98e-13 

1  44e-04  2.49e-l3 

2  88e*04  2.5le-13 
5  76e*04  2  55e-13 

1  15e  +  05  2-11e-13 

2  30e  +  05  2  43e-I3 
4,6le  +  05  188e-13 


lO^  10''  10- 

Averaging  Time,  x,  Seconds 


11.3.  Phase  Noise  Identification 


Consider  the  problem  of  identifying  the  dominant  type  of  phase  noise  in  a  frequency  distribution  amplifier.  Assume 
that  time  domain  stability  measurements  have  been  made  comparing  the  amplifier's  input  and  output  signals.  How 
should  these  data  be  analyzed  to  best  determine  the  noise  type?  Simulated  white  and  flicker  PM  noise  at  equal  levels 
are  analyzed  below  in  several  ways  to  demonstrate  their  ability  to  perform  this  noise  identification. 


Table  27.  Phase  noise  identification 


White  PM 


Flicker  PM 


Examination  of  the  phase  data  is  a  good  first  step.  An  experienced  analyst 
will  immediately  notice  the  difference  between  the  white  and  flicker  noise 
plots,  but  will  find  it  harder  to  quantify  them.  


PHASE  DATA 

simulated  Whtte  PM  NoIm 


4500  5000 


PHASE  DATA 

SInutated  Rlckw  PM  Data 


In  contrast,  the  frequency  data  show  little  noise  type  discrimination, 
because  the  differencing  process  whitens  both.  Examination  of  the 
frequency  data  would  be  appropriate  to  distinguish  between  white  and 
flicker  FM  noise,  however. 


FREQUENCY  DATA 

SImuiatBd  Whne  PM  Ndaa 


FREQUENCY  DATA 

SImulalod  Flicker  PM  Noise 
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The  Allan  deviation  is  not  able  to  distinguish  between  white  and  flicker 
PM  noise.  Both  have  slopes  of  about  1,  as  shown  by  the  superimposed 
noise  fit  lines. 


FREQUENCY  STABILITY 

 SlmJated  WWto  PM  Holaa 


1  i.OOe-09 

2  S.O'c-IO 
4  2  53e-i0 
8  127e-10 

16  6  26e-n 

32  3  07e-n 


;38  7  a3e-r2 

256  3  97e-i2 

512  I95e-12 

102'S  9S3e-ii 


Avera^^  Time.  T.^eoonda 


FREQUENCY  STABILITY 

  Simulated  Fllcnor  PM  Notw 


256  a  20C-12 


lO*  10^  10^ 

Averaging  Time.  T.^Seconds 


The  modified  Allan  deviation,  by  virtue  of  its  additional  phase  averaging, 
is  able  to  distinguish  between  white  and  flicker  PM  noise,  for  which  the 
slopes  are  -1 .5  and  -1 .0,  respectively  


FREQUENCY  STABILITY 

SimUUM  Whn«  PM  NoOa 


AveragLng_  Time,  T.^econds 


FREQUENCY  STABILITY 

^  S^ulaied  Flicker  PM  Nona  


Avera^g  Time,  T,  Seconds 


Even  better  discrimination  is  possible  with  the  autocorrelation  function. 
The  lag  1  ACF  is  0.006  and  0.780  for  these  white  and  flicker  PM  noise 
data,  and  is  able  to  quantitatively  estimate  the  power  law  noise  exponents 
as  +1.99  and  +0.93,  respectively.  That  ID  method  is  quite  effective  even 
for  mixed  noise  types.  


AUTOCORRELATION 

SkrUitM  WHU  PM  Noln 


250      SOO  7S0 


AUTOCORRELATION 

SknuialM  Ricke<  PM  Nolae 


2000    2250  2500 


11.4.  Detection  of  Periodic  Components 

Frequency  stability  data  can  contain  periodic  fluctuations  due  to  external  interference,  environmental  effects,  or  the 
source  itself,  and  it  can  therefore  be  necessary  to  detect  discrete  spectral  components  when  performing  a  stability 
analysis.  This  example  uses  a  set  of  50  000  points  of  t  =  1  second  simulated  white  FM  noise  at  a  level  of  1  x10'"t"'^" 
that  contains  sinusoidal  interference  with  a  period  of  500  s  at  a  peak  level  of  1  x  10"'".  The  latter  simulates 
interference  that  might  occur  as  the  result  of  air  conditioner  cycling.  Several  analytical  methods  are  used  to  detect  this 
periodic  component. 
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Table  28.  Detection  of  periodic  components. 


The  interference  level  is 
too  low  to  be  visible  in  the 
full  frequency  data  plot. 


FREQUENCY  DATA 

Simulated  W  FM  Noise  with  Periodic  Interference 


5.0      10.0      15.0      20.0     25.0     30.0      35.0     40.0     45.0  50.0 

Data  Point  xlOOO 


FREQUENCY  DATA 

Simulated  W  FM  Noise  witti  Periodic  Interference 


By  zooming  in,  one  can  see 
just  a  hint  of  the 
interference  (10  cycles 
over  5000  data  points). 


500      1000     1500     2000    2500     3000     3500     4000    4500  5000 

Data  Point 


FREQUENCY  STABILITY 

Simulated  W  FM  Noise  with  Periodic  Interference 


The  interference  is  quite 
visible  in  an  "all  tau" 
stability  plot  as  a  null  that 
occurs  first  at  an  averaging 
time  of  500  s  (the  period  of 
the  interference).  Here  the 
stability  is  equal  to  the 
underlying  white  FM  noise 
level. 


ISiqfno  for  All  "ajf 

10'  10^  10^  10'' 

Averaging  Tiine^  T,^  Seconds 


10^ 
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The  interference  is  clearly 
visible  in  the  power 
spectral  density  (PSD), 
which  has  a  bright 
component  at  2  mHz, 
correspond-ing  to  the  500  s 
period  of  the  interference. 


POWER  SPECTRUM 

Simulated  W  FM  Noise  with  Periodic  interference 


10-^  10"^  I0-'  10-' 

Fourier  Frequency,  f,  Hz 


AUTOCORRELATION 

Simulated  W  FM  Noise  witli  Periodic  interference 


The  interference  is  less 
visible  in  an 

autocorrelation  plot,  where 
cyclic  behavior  is  barely 
noticeable  with  a  500-lag 
period.  The  PSD  has 
equivalent  information,  but 
its  log  scale  makes  low- 
level  components  more 
apparent. 


250      500      750      1000     1250     1500     1750     2000    2250  2500 

Lag 


It  is  a  good  analysis  policy  to  examine  the  power  spectral  density  when  periodic  fluctuations  are  visible  on  a  stability 
plot  and  periodic  interference  is  suspected. 

11.5.  Stability  Under  Vibrational  Modulation 


This  plot  shows  the  stability 
of  an  oscillator  with  a 
combination  of  white  PM 
noise  and  a  sinusoidal  com- 
ponent simulating 
vibrational  modula-tion. 
Nulls  in  the  Allan  deviation 
occur  at  averaging  times 
equal  to  the  multiples  of  the 
20  Hz  sinusoidal  mod- 
ulation period,  where  the 
stability  is  determined  by 
the  white  PM  noise.  Peaks 
in  the  Allan  deviation  occur 
at  the  modulation  half 
cycles,  and  have  a  t"' 
envelope  set  by  the 
vibrational  phase 
modulation. 


FREQUENCY  STABILITY 

Vibration-Induced  Alian  Deviation  Degradation 


10-'       2  10-' 

Averaging  Time,  t.  Seconds 


103 


1 1 .6.  White  FM  Noise  of  a  Frequency  Spike 


The  Allan  deviation  of  a 
frequency  record  having  a 
large  spike  (a  phase  step)  has 
a  t"""  characteristic  [1].  Thus 
adding  a  single  large  (say 
lO'')  central  outlier  to  the 
1000-point  test  suite  of 
section  12.4  will  give  a  data 
set  with  ay{x)  = 
[(10^)-2/(  1000-1)]"-  = 
3.16386e4,  as  shown  in  this 
stability  plot. 


FREQUENCY  STABILITY 

W  FM  Nol£8  of  Frequency  Spike 


10'        2  10= 
Averaging  Time,  x,  Seconds 


Reference:  C.A.  Greenhall,  "Frequency  Stability  Review,"  Telecommunications 
and  Data  Acquisition  Progress  Report  42-88,  Jet  Propulsion  Laboratory,  Pasadena, 
CA,  Feb.  1987.  

11.7.  Composite  Aging  Plots 


A  composite  plot  showing  the  aging  of  a  population  of  frequency  sources  can  be  an  effective  way  to  present  this 
information  visually,  providing  a  quick  comparison  of  their  behaviors.  The  following  figure  shows  the  stabilization 
period  of  a  production  lot  of  rubidium  clocks. 

Plots  for  40  units  are  shown  with  their  serial  numbers,  all  plotted  with  the  same  scales.  In  particular,  these  plots  all 
have  a  full  x-scale  time  range  of  nine  weeks  (1  week/div)  and  a  full  y-scale  frequency  range  of  1.5  X  10""  (1  X  10' 
'-/div).  A  diagonal  slope  downward  across  the  plot  corresponds  to  an  aging  of  about  -2.4  x  10''Vday.  All  the  data 
have  T=  900  seconds  (15  min).  A  figure  like  this  immediately  shows  that  (a)  all  the  units  have  negative  frequency 
aging  of  about  the  same  magnitude  that  stabilizes  in  about  the  same  way,  (b)  there  are  occasional  gaps  in  some  of  the 
records,  (c)  all  of  the  units  have  about  the  same  short-term  noise,  but  some  of  the  records  are  quieter  than  others  in  the 
longer  term,  (d)  some  of  the  units  take  longer  than  others  to  reach  a  certain  aging  slope. 

The  plots,  although  small,  still  contain  enough  detail  to  allow  subtle  comparisons  very  quickly,  far  better  than  a  set  of 
numbers  or  bar  graphs  would  do.  The  eye  can  easily  see  the  similarities  and  differences,  and  can  immediately  select 
units  based  on  some  criterion,  which  would  be  harder  to  do  using  a  set  of  larger  plots  on  separate  pages.  Closer 
inspection  of  even  these  small  plots  can  reveal  a  lot  of  quantitative  information  if  one  knows  the  scale  factors.  Color 
coding,  although  not  used  here,  could  be  used  to  provide  additional  information. 

These  plots  are  inspired  by  Edward  Tufte's  book  The  Visual  Display  of  Quantitative  Information,  ISBN  978- 
0961392109,  Graphic  Press,  1983. 


104 


105 


12  Software 


Software  is  necessary  to  perform  a  frequency  stability  analysis,  because  those  Specialized  software  is 

calculations  generally  involve  complex  specialized  algorithms  and  large  data  sets  needed  to  perform  a 

needed  to  interactively  perform  and  document  a  complete  stability  analysis.   It  is  frequency  stability  analysis. 
convenient  to  use  an  integrated  software  package  that  combines  all  of  the  required 

analytical  tools,  operates  on  current  computer  hardware  and  operating  systems,  includes  the  latest  analytical 
techniques,  and  has  been  validated  to  produce  the  correct  results. 


12.1.  Software  Validation 


Considerable  effort  is  needed  to  ensure  that  the  results  obtained  from  frequency  stability  analysis  software  are  correct. 
Several  suggested  validation  methods  are  shown.  Mature  commercially  available  software  should  be  used  whenever 
possible  instead  of  developing  custom  software.  User  feedback  and  peer  review  are  important.  There  is  a  continuing 
need  to  validate  the  custom  software  used  to  analyze  time  domain  frequency  stability,  and  the  methods  listed  below 
can  help  ensure  that  correct  results  are  obtained. 

Several  methods  are  available  to  validate  frequency  stability  analysis  software: 

1.  Manual  Analysis:  The  results  obtained  by  manual  analysis  of  small  data  sets  (such  as  in  NBS  Monograph  140 
Annex  8.E)  can  be  compared  with  the  new  program  output.  This  is  always  good  to  do  to  get  a  "feef  for  the 
process. 

2.  Published  Results:  The  results  of  a  published  analysis  can  be  compared  with  the  new  program  output.  One 
important  validation  method  is  comparison  of  the  program  results  against  a  test  suite  such  as  the  one  in 
References  [1]  and  [2].  Copies  of  those  test  data  are  available  on-line  [3]. 

3.  Other  Programs:  The  results  obtained  from  other  specialized  stability  analysis  programs  (such  as  that  from  a 
previous  generation  computer  or  operating  system)  can  be  compared  with  the  new  program  output. 

4.  General  Purpose  Programs:  The  results  obtained  from  industry  standard,  general  purpose  mathematical  and 
spreadsheet  programs  such  as  MathCAD,  Matlab,  and  Excel  can  be  compared  with  the  new  program  output. 

5.  Consistency  Checks:  The  new  program  should  be  verified  for  internal  consistency,  such  as  producing  the  same 
stability  results  from  phase  and  frequency  data.  The  standard  and  normal  Allan  variances  should  be 
approximately  equal  for  white  FM  noise.  The  nonnal  and  modified  Allan  variances  should  be  identical  for  an 
averaging  factor  of  1.  For  other  averaging  factors,  the  modified  Allan  variance  should  be  approximately  one-half 
the  normal  Allan  variance  for  white  FM  noise,  and  the  normal  and  overlapping  Allan  variances  should  be 
approximately  equal.  The  overlapping  method  provides  better  confidence  of  the  stability  estimates.  The  various 
methods  of  drift  removal  should  yield  similar  results. 

6.  Simulated  Data:  Simulated  clock  data  can  also  serve  as  a  useful  cross  check.  Known  values  of  frequency  offset 
and  drift  can  be  inserted,  analyzed,  and  removed.  Known  values  of  power-law  noise  can  be  generated,  analyzed, 
plotted,  and  modeled. 


12.2.  Test  Suites 


Tables  29  and  30  summarize  the  values  for  several  common  frequency  stability  measures  for  both  the  classic  NBS 
data  set  and  a  1000-point  portable  test  suite. 

12.3.  NBS  Data  Set 

A  "classic"  suite  of  frequency  stability  test  data  is  the  set  of  nine  3-digit  numbers  from  Annex  8.E  of  NBS  Monograph 
140  shown  in  Table  29.  Those  numbers  were  used  as  an  early  example  of  an  Allan  variance  calculation.  These 
frequency  data  is  also  normalized  to  zero  mean  by  subtracting  the  average  value,  and  then  integrated  to  obtain  phase 
values.  A  listing  of  the  properties  of  this  data  set  is  shown  in  Table  II.  While  nine  data  points  are  not  sufficient  to 
calculate  large  frequency  averages,  they  are,  nevertheless,  a  very  useful  starting  point  to  verify  frequency  stability 
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calculations,  since  a  small  data  set  can  easily  be  entered  and  analyzed  manually.  A  small  data  set  is  also  an  advantage 
for  finding  "'off-by-one"  errors  where  an  array  index  or  some  other  integer  variable  is  slightly  wrong  and  hard  to 
detect  in  a  larger  data  set. 


Table  29.  NBS  Monograph  140,  Annex  8.E  Test  Data 


# 

Frequency 

Normalized  frequency 

Phase  (t=1) 

1 

892 

103.1 1 1 1 

0.00000 

2 

809 

20.11111 

103.1 1 1 1 1 

823 

34.11111 

123  22222 

4 

798 

9.1 11 11 

157.33333 

5 

671 

-1 17.88889 

166.44444 

6 

644 

-144.88889 

48.55555 

7 

883 

94.11111 

-96.33333 

8 

903 

114.11111 

-2.22222 

9 

677 

-1 1 1.88889 

111.88889 

10 

0.0000 

Table  30.  NBS  Monograph  140,  Annex  8.E  Test  Data  Statistics 


Averaging  Factor 

1 

2 

#  Data  Points 

9 

4 

Maximum 

903 

893.0 

Minimum 

644 

657.5 

Average 

788.8889 

802.875 

Median 

809 

830.5 

Linear  Slope 

-10.20000 

-2.55 

Intercept 

839.8889 

809.25 

Standard  Deviation^^' 

100.9770 

102.6039 

Normal  Allan  Deviation 

91.22945 

115.8082 

Overlapping  Allan  Dev 

91.22945 

85.95287 

Modified  Allan  Dev 

91.22945 

74.78849 

Time  Deviation 

52.67135 

86.35831 

Hadamard  Deviation 

70.80607 

1 16.7980 

Overlap  Hadamard  Dev 

70.80607 

85.61487 

Hadamard  Total  Dev 

70.80607 

91.16396 

Total  Deviation 

91.22945 

93.90379 

Modified  Total  Dev 

75.50203 

75.83606 

Time  Total  Deviation 

43.59112 

87.56794 

Note:  [a]  Sample  (not  population)  standard  deviation. 

12.4.  1000-Polnt  Test  Suite 

The  larger  frequency  data  test  suite  used  here  consists  of  1000  pseudo-random  frequency  data  points.  It  is  produced 
by  the  following  prime  modulus  linear  congruential  random  number  generator: 

A7,^,  =16807^,  Mod  2147483647.  (73) 


107 


This  expression  produces  a  series  of  pseudo-random  integers  ranging  in  value  from  0  to  2  147  483  646  (the  prime 
modulus,  2^'-l,  avoids  a  collapse  to  zero).  When  started  with  the  seed  no  =  1  234  567  890,  it  produces  the  sequence  ni 
=  395  529  916,  n.  =  1  209  410  747,  nj  =  633  705  974,  etc.  These  numbers  may  be  divided  by  2  147  483  647  to  obtain 
a  set  of  normalized  floating-point  test  data  ranging  from  0  to  1.  Thus  the  normalized  value  of  no  is  0.5748904732.  A 
spreadsheet  program  is  a  convenient  and  reasonably  universal  way  to  generate  these  data.  The  frequency  data  may  be 
converted  to  phase  data  by  assuming  an  averaging  time  of  1,  yielding  a  set  of  1001  phase  data  points.  Similarly, 
frequency  offset  and/or  drift  terms  may  be  added  to  the  data.  These  conversions  can  also  be  done  by  a  spreadsheet 
program. 

The  values  of  this  data  set  will  be  uniformly  distributed  between  0  and  1.  While  a  data  set  with  a  normal  (Gaussian) 
distribution  would  be  more  realistic,  and  could  be  produced  by  summing  a  number  of  independent  uniformly 
distributed  data  sets,  or  by  the  Box-MuUer  method  [5],  this  simpler  data  set  is  adequate  for  software  validation. 


Table  31.  1000-point  frequency  data  set. 


Averaging  factor 

1 

10 

100 

#  Data  Points 

1000 

100 

10 

Maximum 

9.957453e-01 

7.003371e-01 

5.489368e-01 

Minimum 

1.371760e-03 

2.545924e-01 

4.533354e-01 

Average'^' 

4.897745e-01 

4.897745e-01 

4.897745e-01 

Median 

4.798849e-01 

5.047888e-01 

4.807261e-01 

Linear  Slope''^'^' 

6.4909  lOe-06 

5.979804e-05 

1.056376e-03 

 Fl   

Intercept^'"' 

4.865258e-01 

4.867547e-01 

4.839644e-01 

Bisection  Slope'^' 

-6.104214e-06 

-6.104214e-05 

-6.104214e-04 

1st  Diff  Slope'^' 

1.517561e-04 

9.648320e-04 

1.01 1791e-03 

Log  Fitt^''',  a= 

5.577220e-03 

7  1 38988P-03 

y(t)=a-ln  (bt+l)+c,  b= 

9.737500e-01 

4.594973e+00 

1.420429e+02 

y'(t)=ab/(bt+l),  c= 

4.570469e-01 

4.631 172e-01 

4.442759e-01 

Slope  at  end 

5.571498e-06 

5.237080e-05 

7.133666e-04 

Standard  Dev'" 

2.884664e-01 

9.296352e-02 

3.206656e-02 

Normal  Allan  Dev'*^ 

2.9223 19e-01 

9.965736e-02 

3.897804e-02 

Overlap  Allan  Devf*"' 

2.9223 19e-01 

9.159953e-02 

3.241343e-02 

Mod  Allan  Dev'='^' 

2.9223  19e-01 

6.172376e-02 

2.1 7092  le-02 

Time  Deviation'*^' 

1.687202e-01 

3.563623e-01 

1.253382e-00 

Hadamard  Deviation 

2.943883e-01 

1.052754e-01 

3.910861e-02 

Overlap  Had  Dev 

2.943883e-01 

9.581083e-02 

3.237638e-02 

Hadamard  Total  Dev 

2.943883e-01 

9.614787e-02 

3.058103e-02 

Total  Deviation 

2.9223 19e-01 

9.134743e-02 

3.406530e-02 

Modified  Total  Dev 

2.418528e-01 

6.4991 61e-02 

2.287774e-02 

Time  Total  Deviation 

1.396338e-01 

3.752293e-01 

1.320847e-00 

Notes: 

[a]  Expected  value  =  0.5. 

[b]  All  slopes  are  per  interval. 

[c]  Least  squares  linear  fit. 

[d]  Exact  results  will  depend  on  iterative  alg 

orithm  used.  Data  not  suited  to  log  fit. 

[e]  Sample  (not  population)  standard  deviation.  Expected 

value  =   1/V12  = 

0.2886751. 

[f]  Expected  value  equal  to  standard  deviation  for  white  FM  noise. 

[g]  Equal  to  normal  Allan  deviation  for  averaging  factor  =  1 . 

[h]  Calculated  with  listed  averaging  factors  from  the  basic  i  =  1 

data  set. 
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Table  32.  Error  bars  for  n=1000  point  =1  data  set  with  avg.  factor=10. 


Allan  Dev  Type 
#Pts 


Sigma  Value 
Noise  Type  &  Ratio 


Conf  Interyal 
Remarks 


95%  CF 


Normal 
99 


9.965736e-02 

W  FM^'l  Bl=0.870 


CI=  8.7138706-03^^'^ 
ilaClt^. 


Overlapping 
981 


9.159953e-02 
W  FM 

#     df-  146.177 


Max  ay(T)=  1.01 4923 e-Ol'^^l 
Max  Gy(T)=  1. 03520 le-Ol''^' 
Min  gJx)=  8.2239426-02'^' 


19.07 
14.45 
81.34 


Modified^'^l 
972 


6.1723766-02 

W  FM'"',  R(n)=0.384 

#x'  df=  94.620 


Max  ay(T)=7.04441 26-02'*^' 
Max  a,(T)=  7.2249446-02''^' 
Min  av(i)-  5.4199616-02''^' 


72.64 
69.06 
122.71 


Notes: 

[a]  Theoretical  B  1  =  1. 000  for  W  FM  noise  and  0.667  for  F  and  W  PM  noise. 

[b]  Simple,  noise-independent  CI  estimate  =ay(T)/VN=l  .001 594e-02. 

[c]  This  CI  includes  K(a)  factor  that  depends  on  noise  type: 


Noise 

a 

K(a) 

W  PM 

2 

0.99 

FPM 

1 

0.99 

W  FM 

0 

0.87 

F  FM 

-1 

0.77 

RW  FM 

.2 

0.75 

[d]  BW  factor  2T:fhTo=  10.  Applies  only  to  F  PM  noise. 

[e]  Theoretical  R(n)  for  W  FM  noise  =  0.500  and  0.262  for  F  PM  noise. 

[f]  Double-sided  68  %  confidence  interval:  p  =  0.158  and  0.842. 

[g]  Single-sided  95  %  confidence  interval:  p  =  0.950. 

[h]  Double-sided  95  %  confidence  interval:  p  =  0.025  and  0.975. 


12.5.  IEEE  Standard  1139-1999 

IEEE  Standard  1139-1999,  IEEE  Standard  Definitions  of  Physical  Quantities  for  Fundamental  Frequency  and  Time 
Metrology  -  Random  Instabilities  contains  several  examples  of  stability  calculations.  Annex  C.2  contains  an  example 
of  an  Allan  deviation  calculation.  Annex  C.3  has  an  example  of  a  modified  Allan  deviation  calculation.  Annex  C.4 
has  an  example  of  a  total  deviation  calculation,  and  Annex  D  contains  examples  of  confidence  interval  calculations. 
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13  Glossary 


The  following  terms  are  used  in  the  field  of  frequency  stability  analysis 


Aging 

Allan  variance 
Averaging 
Averaging  time 
BW 

Confidence  limit 
Drift 

Frequency  data 
Hadamard  variance 

£(f) 
MJD 

Modified  sigma 
MTIE 

Normalize 
Phase  data 

Phase  noise 

Sampling  time 

Sigma 

Slope 

SSB 

Sf(f) 


The  change  in  frequency  with  time  due  to  internal  effects  within  the  device. 

The  two-sample  variance  CT^y(T)  commonly  used  to  measure  frequency  stability. 

The  process  of  combining  phase  or  frequency  data  into  samples  at  a  longer  averaging  time. 

See  Tau. 

Bandwidth,  hertz. 

The  uncertainty  associated  with  a  measurement.  Often  a  68  %  confidence  level  or  error 
bar. 

The  change  in  frequency  with  time  due  to  all  effects  (including  aging  and  environmental 
sensitivity). 

A  set  of  fractional  frequency  values,  y[i],  where  i  denotes  equally-spaced  time  samples. 

A  three-sample  variance,  HVAR,  that  is  similar  to  the  two-sample  Allan  variance.  It  uses 
the  second  differences  of  the  fractional  frequencies,  and  is  unaffected  by  linear  frequency 
drift. 

£(f)  =  101og['/2  •  Si,{f)],  the  ratio  of  the  SSB  phase  noise  power  in  a  1  Hz  BW  to  the  total 
carrier  power,  dBc/Hz.  Valid  when  noise  power  is  much  smaller  than  the  carrier  power. 

The  Modified  Julian  Date  is  based  on  the  astronomical  Julian  Date,  the  number  of  days 
since  noon  on  January  1,  47 13  BC.  The  MJD  is  the  Julian  Date  -  2  4000  000.5. 

A  modified  version  of  the  Allan  or  total  variance  that  uses  phase  averaging  to  distinguish 
between  white  and  flicker  PM  noise  processes. 

The  maximum  time  interval  error  of  a  clock. 

To  remove  the  average  value  from  phase  or  frequency  data. 

A  set  of  time  deviates,  x[i]  with  units  of  seconds,  where  i  denotes  equally-spaced  time 
samples.  Called  "phase"  data  to  distinguish  them  from  the  independent  time  variable. 

The  spectral  density  of  the  phase  deviations. 

See  Tau. 

The  square  root  or  deviation  of  a  variance,  often  the  two-sample  or  Allan  deviation,  CTy(i). 
The  change  in  frequency  per  tau  interval. 
Single  sideband. 

The  one-sided  spectral  density  of  the  phase  deviations,  radVHz. 
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Sx(f) 
Sy(f) 

Tau 
Total 

Theo1 

TheoBR 
TheoH 

TIE 

Total  variance 
x(t) 

y(t) 


The  one-sided  spectral  density  of  the  time  deviations,  secVHz. 

The  one-sided  spectral  density  of  the  fractional  frequency  deviations,  1/Hz. 

The  interval  between  phase  measurements  or  the  averaging  time  used  for  a  frequency 
measurement. 

A  variance  using  an  extended  data  set  that  provides  better  confidence  at  long  averaging 
times. 

Theoretical  variance  #1,  a  variance  providing  stability  data  out  to  75  %  of  the  record 
length. 

A  biased-removed  version  of  Theol. 

A  hybrid  combination  of  TheoBR  and  the  Allan  variance. 

The  time  interval  error  of  a  clock.  Can  be  expressed  as  the  rms  time  interval  error  TIE  rms 
or  the  maximum  time  interval  error  MTIE. 

A  two-sample  variance  similar  to  the  Allan  variance  with  improved  confidence  at  large 
averaging  factors. 

The  instantaneous  time  deviation  from  a  nominal  time,  x(t)  =  <!f>{t)/2nvo,  s,  where  v,,  is  the 
nominal  frequency,  hertz.  This  dependent  time  variable  is  often  called  "phase''  to 
distinguish  it  from  the  independent  time  variable  t. 

The  instantaneous  fractional  frequency  deviation  from  a  nominal  frequency,  y(t)  = 
[v(t)-vo]/vo]  =  x'(t),  where  vq  is  the  nominal  frequency. 
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Random  Instabilities,"  IEEE  Std  1 139-1999,  July  1999. 
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Communication  -  International  Telecommunications  (CCITT)  Union,  Geneva,  Switzerland. 
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